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ABSTRACT 

Sr,  Nd  and  Pb  isotopic  analysesc^f  477  samples  representing  30  islands  or 
island  groups,  3  seamounts  or  seamount  chains,  2  oceanic  ridges  and  1  oceanic 
plateau  (for  a  total  of  36  geographic  features]  are  compiled  to  form  a 
comprehensive  oceanic  island  basalt  fOlB]  data  set.  These  samples  are 
supplemented  by  90  selected  midt50cean  ridge  basalt  [MORE]  samples  to  give 
adequate  representation  to  MORE  as  an  oceanic  basalt  end-member.  This 
comprehensive  data  set  is  used  to  infer  information  about  the  Earth's  mantle. 
Principal  component  analysis  of  the  OIB-f-MORB  data  set  shows  that  the  first 
three  principal  components  account  for  97.5%  of  the  variance  of  the  data.  Thus, 
only  four  mantle  end-member  components  (EMI,  EMII,  HIMU  and  DMM]  are 
required  to  completely  encompass  the  range  of  known  isotopic  values.  Each 
sample  is  expressed  in  terms  of  percentages  of  the  four  mantle  components, 
assuming  linear  mixing. ,  There  is  significant  correlation  between  location  and 
isotopic  signature  within  geographic  features,  but  not  between  them,  so 
discrimination  analysis  of  the  viability  of  separating  the  oceanic  islands  into  those 
lying  inside  and  outside  Hart's  (1984,  1988)  DUPAL  belt  is  performed  on  the 
feature  level  and  yields  positive  results. 

/ 

A  'continuous  layer  modeF  is  applied  to  the  mantle  component  percentage 
data  to  solve  for  the  spherical  harmonic  coefficients  using  approximation 
methods.  Only  the  degrees  0-5  coefficiwits  can  be  solved  for  since  there  are  only 
36  features.  The  EMI  and  HIMU  percentage  data  sets  must  be  filtered  to  avoid 
aliasing.  D^to  the  nature  of  the  data,  the  coefficients  must  be  solved  for  using 
singular  valu^^lecomposition  (SVD),  versus  the  least  squares  method.  The  F-test 
provides  an  objective  way  to  estimate  the  number  of  singular  values  to  retain 
when  solving  with  SVD.  With  r'*‘^pe'?t  to  the  behavior  of  geophysics  control  data 
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sets,  only  the  degree  2  spherical  harmonic  coefficients  for  the  mantle  components 
can  be  estimated  with  a  reasonable  level  of  confidence  with  this  method. 

Applying  a  "delta-function  model"  removes  the  problem  of  aliasing  and 
simplifies  the  spherical  harmonic  coefficient  solutions  from  integration  on  the 
globe  to  summation  over  the  geographic  features  due  to  the  propeities  of  delta- 
functions.  With  respect  to  the  behavior  of  geophysics  control  data  sets,  at  least 
the  degree  2  spherical  harmonic  coefficients  for  the  mantle  components  can  be 
estimated  with  confidence,  if  not  the  degrees  3  and  4  as  well.  Delta-function 
model  solutions  are,  to  some  extent,  controlled  by  the  nonuniform  feature 
distribution,  while  the  continuous  layer  model  solutions  are  not. 

The  mantle  component  amplitude  spectra,  for  both  models,  show  power  at 
all  degrees,  with  no  one  degree  dominating.  The  DUPAL  components  [EMI, 
EMU  and  HIMU],  for  both  models,  correlate  well  with  the  degree  2  geoid, 
indicating  a  deep  origin  for  the  components  since  the  degrees  2-3  geoid  is 
inferred  to  result  from  topography  at  the  core-mantle  boundary.  The  DUPAL 
and  DMM  components,  for  both  models,  correlate  well  [and  negatively]  at  degree 
3  with  the  velocity  anomalies  of  the  Clayton-Comer  seismic  tomography  model 
in  the  2500-2900  km  depth  range  [immediately  above  the  core-mantle  boundary]. 
The  EMll  component  correlates  well  [and  positively]  at  degree  5  with  the 
velocity  anomalies  of  the  Clayton-Comer  model  in  the  700-1200  km  depth  range, 
indicating  a  subduction  related  origin.  Similar  positive  correlations  for  the  geoid 
in  the  upper  lower  mantle  indicate  that  subducted  slabs  extend  beyond  the  670 
km  seismic  discontinuity  and  support  a  whole-mantle  convection  model. 

Thesis  Supervisor;  Dr.  Stanley  R.  Hart 

Title:  Senior  Scientist,  Woods  Hole  Oceanographic  Institution 
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CHAPTER  1 
INTRODUCTION 

PREVIOUS  WORK 

That  the  Earth's  mantle  is  heterogeneous  is  no  longer  a  subject  of 
controversy  among  geochemists,  but  the  composition,  the  location  and  the 
geometry  of  these  heterogeneities  is  very  much  in  question.  Direct  sampling  is 
not  an  option  for  studying  the  chemistry  of  most  of  the  mantle,  so  products  of 
indirect  sampling,  such  as  oceanic  island  basalts  fOIB'sj  and  mid-ocean  ridge 
basalts  [MORE],  are  invaluable  for  revealing  the  nature  of  the  inaccessible 
mantle.  Though  the  OIB's  may  be  contaminated  by  interactions  with  the 
lithosphere  or  may  sample  large  vertical  sections  of  the  mantle,  they  still  retain 
the  signature  of  their  original  source. 

Using  various  statistical  methods  and  models,  previous  workers  have 
defined  what  they  believe  to  be  the  number  of  mantle  component  end-members 
required  to  represent  the  variation  in  the  oceanic  mantle  data  lOlB-t-MORB). 
Early  on,  Zindler  et  al.  (1982)  used  factor  analysis  to  evaluate  the  oceanic  data  in 
five  dimensions.  Their  analysis  indicated  that  the  oceanic  data  define  a  plane  (the 
"mantle  plane"],  described  by  the  mixing  of  three  chemically  independent 
components,  two  undifferentiated  or  slightly  enriched  mantle  components  and 
one  MORB-type  or  depleted  mantle  component. 

Other  workers  have  chosen  five  groups  or  components  to  represent  the 
data.  Using  a  series  of  two-dimensional  isotopic  plots.  White  (1985)  divided  the 
oceanic  data  into  five  distinct  basalt  groups  [MORB,  St.  Helena,  Kerguelen, 
Society,  and  Hawaii].  He  concedes  that  the  five  groups  may  be  end-members 
which  mix  to  form  intermediate  compositions,  but  he  believes  that  each  group 
either  represents  a  distinct,  internally  homogeneous  reservoir  or  that  each  group 


is  composed  of  a  number  of  isotopically  similar  reservoirs.  Likewise,  Li  e!  al. 
(1991)  proposed  fives  extremes,  using  non-linear  mapping:  Atlantic  MORB 
[DMM],  St.  Helena  [HIMU],  Walvis  [EMI],  Samoa  [EMU]  and  [EMlllj.  Non¬ 
linear  mapping  approximately  preserves  the  geometric  structure  of  the  data  by 
maintaining  interpoint  distances.  Four  of  the  five  extremes  of  Li  et  al.  (1991 ) 
are  based  solidly  on  samples  trends  from  islands,  but  the  D5  extreme  is  based 
only  on  that  one  sample.  More  data  is  needed  to  substantiate  their  fifth  extreme. 

By  far  the  majority  of  analyses  indicate  the  existence  of  four  end-member 
components  for  the  oceanic  mantle  data.  Using  two-dimensional  plots,  Zindler 
and  Hart  (1986)  defined  the  following  four  end-member  components:  depleted 
MORB  mantle  [DMM],  high  U/Pb  mantle  [HIMU],  and  two  enriched  mantle 
components  [EMI  and  EMU],  with  possibly  two  other  components  prevalent 
mantle  composition  [PREMA]  and  bulk  silicate  Earth  [BSE],  Eigenvector 
analyses  by  Allegre  et  al.  (1987)  agree  with  the  four  component  model  of 
Zindler  and  Hart  (1986).  The  four  extremes  of  Allegre  et  al.  (1987)  are 
[correspond  to]:  extreme  MORB  [DMM];  St.  Helena,  Tubuai  and  Mangai  islands 
[HIMU];  Kerguelen,  Gough,  Tristan  da  Cunha  and  Raratonga  islands  |EM1];  and 
Sao  Miguel  and  Atui  islands  [EMU].  Hart  (1988),  using  an  augmented  data  set 
and  two-dimensional  plots,  concluded  that  the  four  end-members  proposed  by 
Zindler  and  Hart  (1986)  are  valid  representations  of  the  extremes  of  the  oceanic 
data.  He  resolves  White's  (1985)  groupings  into  his  own  four  component  system 
as  follows  [White  =  Hart]:  MORB  =  DMM,  Society  =  EMU,  St.  Helena  =  HIMU. 
Hawaii  =  EMI,  with  the  suggestion  that  White's  fifth  group,  Kerguelen,  is  a 
mixture  of  EMI  and  EMIL  In  addition,  Li  et  al.  (1991)  also  noted  a  tetrahedral 
structure  to  the  data,  when  using  factor  analysis  with  varimax  rotation,  with  the 
following  four  extremes:  Atlantic  MORB  [OMM],  Mangaia  [HIMU],  Samoa 
[EMU]  and  WaIvis[EMI]. 


One  scenario  for  the  genesis  of  the  three  unusual  mantle  components  is  put 
forth  by  Hart  (1988).  He  proposes  that  HIMU,  enriched  in  U,  is  probably 
generated  by  intra-mantle  metasomatism,  that  EMI  corresponds  to  a  slightly 
modified  bulk-earth  compositon  and  that  EMU  can  be  explained  by  the  recycling 
of  sediments  during  subduction.  The  proposed  formation  mechanisms  in  no  way 
limit  the  geometry  of  the  mantle  needed  to  generate  the  heterogeneities  and,  as 
such,  a  wide  variety  of  models  have  been  proposed.  A  whole  mantle  convection 
model  might  portray  the  enriched  mantle  components  as  blobs  floating  around  in 
a  depleted  mantle  matrix  (Zindler  and  Hart,  1986)  or  perhaps  as  an  accutnulated 
layer  of  subducted  oceanic  crust  and  sediment  at  the  core-mantle  boundary  that 
reaches  the  surface  in  mantle  plumes  (Hofmann  and  White.  1982).  A  two-layer 
convection  model  might  rely  on  a  depleted  upper  mantle  feeding  the  mid-ocean 
ridges  and  an  enriched  low'cr  mantle  feeding  oceanic  islands  via  mantle  plumes 
(Dupre  and  Allegre,  1983)  or  require  a  depleted  upper  mantle,  a  primitive  lower 
mantle  and  an  accumulated  layer  of  subducted  oceanic  crust  and  sediment  at  the 
67U  km  discontinuity  that  supplies  the  enriched  components  via  mantle  plumes 
(White,  1985;  Allegre  and  Turcotte,  1985).  Anderson  (1985)  even  proposes  a 
three-layer  convective  model  with  the  geochemical  contrasts  occurring  vS.'.ly  in 
the  upper  mantle  with  a  depleted  lower  part  that  supplies  the  mid-ocean  ridges 
and  an  enriched  upper  part  from  subduction  of  oceanic  crust  and  sediment. 

A  deep  origin  for  the  enriched  components  is  indicated  by  Hart's  (1984) 
large-scale  isotopic  anomaly,  the  DUPAL  anomaly,  characterized  by  the 
concentration  of  the  enriched  mantle  components  in  a  band  from  2°  S  to  60°  S. 
Qualitatively,  countours  of  the  anomaly  criteria  [A7/4,  A8/4  and  ASr  (Hart, 

1984)1  correspond  to  lone-w’avelength  [and  thus  deep]  geophysical  quantities 
(Hart,  1988).  Other  researchers  oppo,se  this  deep  origin  interpretation,  citing  die 
nonuniform  distribution  of  hotspots  as  the  reason  for  the  pattern  (White.  1985) 


or  arguing  that  the  DUF^AL  compositions  occur  in  scattered  locations  and  do  not 
can’cr  a  coherent  geographic  area  (Allegre  et  ai.  1987). 

The  purpose  of  this  thesis  is  three-fold:  (1)  to  address  once  again  the  issue 
of  the  number  of  mantle  end-member  components  needed  to  represent  the 
oceanic  mantle  data,  H)  to  statistically  test  the  viability  of  the  DUF^AI.  distinction 
as  a  means  of  characterizing  the  OIB  data  and  (3)  to  try  to  pinpoint  the  source 
depth  of  the  enriched  mantle  components  by  expanding  their  relative  abundances 
ill  spherical  harmonics  and  comparing  their  expansions  to  those  of  knou  n 
geophysical  quantities. 


DATA 

The  majority  of  this  study  focuses  on  Sr.  Nd  and  Pb  isotopic  analyses  of 
volcanic  rocF:s  from  oceanic  islands,  seamounts,  lidges.  and  plateaus.  All  of 
these  geographic  features  overlie  oceanic  crust,  with  the  exception  of  Nunivak 
Island  on  the  Alaskan  Continental  Shelf,  and  none  of  them  is  directly  associated 
with  seafloor  spreading,  with  the  exception  of  Iceland,  which  has  a  mixture  of 
mid-ocean  ridge  and  hotspot  influences.  Essentially,  the  data  set  is  that  compiled 
by  Zindler  et  al.  (1982)  and  later  augmented  by  Hart  (1988),  with  some 
additional  recent  analyses  (Appendix).  Samples  in  the  data  set  are  mainly  btisali. 
with  some  gabbros  and  trachybasalts;  trachytes  and  other  silica-rich  rocks 
relative  to  basalt  | roughly  SiO-)  >  50% |  are  excluded.  The  majority  of  the 

samples  are  of  Cenozoic  age.  with  the  exception  of  the  Walvis  Ridge,  Rio  Grande 
Rise  and  New  England  Seamounts  samples,  with  ages  up  to  100  Ma.  If  a  choice 
is  given,  analyses  of  leached  samples  are  preferred  over  analyses  of  unleacheH 
samples.  In  addition,  only  single  samples  for  which  there  are  Sr,  Nd  and  Pb 
analyses  are  included.  Fmr  ctinsistency.  Sr  data  is  adjusted  to  0.7()8(K)  |li&A 
sttindardl  or  0.71022  |NBS  SRM  987  standard]  and  Nd  data  is  adjusted  to 


0. 512640  [BCR-1  standard]  or  0.51 1862  [La  Jolla  standard]  or  0.51 1296  ]Spcx 
standard  ]. 

In  this  data  set,  referred  to  as  the  OIB  data  set.  there  are  477  samples 
representing  30  islands  or  island  groups,  3  seamounts  or  seamount  chains,  2 
aseismic  oceanic  ridges  and  1  oceanic  plateau  (Figure  1.1  and  Table  1.1).  J'he 
isotopic  means  and  standard  deviations  for  the  OIB  data  are  listed  in  Table  1 .2. 

Since  MORB  is  considered  to  be  one  of  the  mantle  component  end- 
members  (Zindler  et  aL,  1982;  White,  l9S5;  Zindler  and  Hart,  1986;  Hart. 

1988),  any  attempt  to  choose  end-members  should  include  MORB  data,  b'or  this 
reason,  a  second  data  set  is  created  using  the  OIB  data  and  a  selection  of  90 
.MORB  samples  (Appendix),  the  OIB+MORB  data  set  (Table  1.3).  The  criteria 
for  choosing  OIB  samples  applies  to  the  MORB  samples  as  well.  Isotopic  means 
and  standard  deviations  for  the  OIB-i-MORB  data  are  listed  in  Table  1.2. 

ORGANIZATION 

The  main  thust  of  this  work  is  to  characterize  the  OIB  data  and  to  search 
for  possible  correlations  between  the  geochemical  signatures  of  OlB's  and 
geophysical  quantities,  such  as  the  geoid  and  .seismic  tomography,  that  might  help 
pinpoint  the  depth].sl  of  the  OIB  reservoirjsj. 

Chapter  2  explores  the  nature  of  the  OIB  isotope  data.  With  the  help  of 
principal  component  analysis,  the  data  is  expressed  in  terms  of  percentages  of 
four  mantle  component  end-members.  Spatial  correlation  testing  reveals  the 
relationship  between  geographic  distance  from  island  to  island  and  feature  to 
feature  and  the  "isotopic  distance"  between  samples.  Discrimination  analysis, 
both  nearest-neighbor  and  graphical,  is  u.sed  to  test  the  viability  of  separating  the 
oceanic  islands  into  two  groups,  inside  and  outside  the  DUPAL  belt. 


Chapter  3  applies  a  "continuous  layer  model"  to  the  mantle  component 
data,  as  an  assumed  geometry  for  the  OIB  reservoir,  in  order  to  solve  for  the 
spherical  harmonic  coefficients.  The  problem  of  aliasing  is  addressed  with  the 
relationship  of  variation  in  mantle  components  to  distance  between  features. 
Approximation  methods  are  used  to  solve  for  the  coefficients.  Geophysical  data 
sets  are  constructed,  using  GEM-L2  geoid  coefficients,  to  serve  as  controls 
against  which  to  judge  the  success  of  the  approximation  methods. 

Chapter  4  applies  a  "delta-function  model"  to  the  mantle  component  data 
to  provide  a  mathematically  more  robust  solution  for  the  spherical  harmonic 
coefficients.  The  delta-function  approximation  removes  the  problem  of  aliasing, 
but  generates  a  solution  dependent  upon  feature  location.  The  same  geophysical 
data  sets  are  used  again  to  judge  the  success  of  the  delta-function  approximation. 

Chapter  5  compares  the  mantle  component  spherical  harmonic  solutions 
for  the  two  models  in  terms  of  their  amplitude  spectra,  how  well  they  correlate 
with  the  geoid,  how  they  are  affected  by  the  nonuniform  feature  distribution  and 
how  well  they  correlate  with  the  Clayton-Comer  seismic  tomography  model. 

The  implications  of  these  results  and  recommendations  for  further  research  are 
discussed. 


Table  1.1.  Geographic  features  represented  in  the  OIB  data  set,  with  their 
components,  number  of  samples  (in  braces]  and  references  indicated. 


Feature 

Components 

References  ’ 

Ascension  [5]  7,34,35 

Amsterdam/St.  Paul  [11] 

Amsterdam  [5]  38 

St.  Paul  |6]  38 

Azores  [6] 

Faialfl]  22 

Sao  Miguel  [5]  1,8 

Balleny  [3]  19 

Cameroon  Line  [18] 

Bioko  15]  17,18 

Pagalu  [  1  ]  18 

Principe  (3]  1 8 

Sao  Tome  [9]  17,18 

Cape  Verde  Islands  [41  [ 

Fogo  [6]  14 

Maio[9]  8,14 

SaoAntao[10]  8,14 

SaoTiago[13]  14 

Sao  Vincente  [3]  14 

Christmas  [  13]  19 

Cocos  [3 1  3 

Comores  Archipelago  [14]  38 

Cook- Austral  Islands  [261 

Aitutaki  [4]  1,21 

Atui  [6]  1,21 

Mangia[5[  1,21 

Mauke[3)  1,21 

Raratonga  [8]  1,21,23 

Crozet  Islands  [9]  38 

Fernando  de  Noronha  [16]  1,13 

Galapagos  Islands  [11]  39 

Gough  [2]  1 

Hawaiian  Islands  [73] 

Hawaii  [14]  28,32 

Kahoolawe  [13]  37 

Kauai  [2]  28 

Lanai  [4]  37 


Table  1.1.  Continued. 
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Feature 


Components 


References 


Hawaiian  Islands  [73] 

Loihi  [15]  27 

Maui  [3]  28 

Molokai  [5]  28 

Oahu  [17]  29 

Iceland  [7]  20 

Juan  Fernandez  Islands  [4]  15 

Kerguelen  Plateau  [41 1 

Heard  Island  [9]  2,30 

Kerguelen  Island  [20]  12,30,38 

Kerguelen  Plateau  [12]  26,36 

Louisville  Seamount  Chain  [4]  6 

Marion/Prince  Edward  [4]  19 

Marquesas  Archipelago  [11]  1 0, 1 1 ,33 

Mascareignes  [8] 

Mauritius  [1]  1 

Reunion  [7]  38 

New  England  Seamounts  [6]  31 

■  Nunivak  [2]  25 

Pitcairn  [19]  41 

Ponape[l]  19 

Sala  Y  Gomez  [  1  ]  ? 

Samoa  Islands  [34] 

Manu'a  [4]  42 

Savai'i  [8]  42 

Tutuila  [9]  23,42 

Upolu[13]  23,42 

San  Felix/San  Ambrosio  15| 

San  Felix  [4]  15 

San  Ambrosio  [  1  ]  15 

Shimada  Seamount  [1  ]  16 

Society  Ridge  [9] 

Mehetia  [2]  9 

Moua  Pihaa  [  1  ]  9 

Tahaa[l]  40 

Teahitia  [4]  9 

dredge  [  1  ]  9 

7,1,4,7,22 


St.  Helena  [31] 


Table  1.1.  Continued. 


Feature 

Components 

References 

Trinidade  [1] 

1 

Tristan  de  Cunha  [5] 

7,22 

Tubuai-Austral  Islands  [22] 

Marotiri  [1] 

5 

Raevavae [ 1 ] 

1 

Rapa  [3] 

5,23 

Rimatara  [4] 

7,21 

Rurutu  [4] 

21,23 

Tubuai  [9] 

5 

Walvis  Ridge  [10] 

24 

Mn  the  reference  column,  a  "?"  indicates  a  sample  with  an  unknown  reference. 

Reference  guide:  ( 1  ]  Allegre  et  ai,  1987;  [2]  Barling  and  Goldstein.  1 990;  |3| 
Castillo  etai,  1988;  [4]  Chaffey  et  al.,  1989;  [5]  Chauvel  eiai,  1991;  16]  Cheng 
et  ai,  1988;  [7]  Cohen  and  O'Nions,  1982a;  [8]  Davies  et  ai,  1989;  [91  Devey  et 
ai,  1990;  [10]  Duncan  et  al.,  1986;  [11]  Dupuy  etai,  1987;  [12]  Gautier  c/ a/., 
1990;  [13]  Gerlach  et  ai,  1987;  [14]  Gerlach  etai,  1988;  [15]  Gerlach  et  ai, 
1986;  [16]  Graham,  1987;  [17]  Halliday  ef  a/.,  1990;  [18]  Halliday  etai,  1988; 
[19]  Hart,  1988;  [20]  Hart,  unpublished;  [21]  Nakamura  and  Tatsumoio,  1988; 

[22]  Newsom  et  al.,  1986;  [23]  Palacz  and  Saunders,  1986;  [24]  Richardson  et  ai, 
1982;  [25]  Roden,  1982;  [26]  Salters,  1989;  [27]  Staudigel  etai.,  1984;  [28]  Stille 
et  ai,  1986;  [29]  Stille  et  ai,  1983;  [30]  Storey  et  ai,  1988;  [31]  Taras  and  Hart, 
1987;  [32]  Tatsumoto,  1978;  [33]  Vidal  etai,  1984;  [34]  Weis,  1983;  [35]  Weis  et 
ai,  1987;  [36]  Weis  et  ai,  1989;  [37]  West  et  ai,  1987;  [38]  White,  unpublished; 
[39]  White  and  Hofmann,  1982;  [40]  White  et  ai,  1989;  [41  ]  Woodhead  and 
McColloch,  1989;  [42]  Wright  and  White,  1987. 
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Table  1.2.  Isotopic  means  and  standard  deviations^  for  the  OIB  and  the 
OIB+MORB  data  sets. 


Sr 

Nd 

6/4Pb 

7/4  Pb 

8/4Pb 

0IB2 

Mean 

0.703943 

0.512825 

19.065 

15.586 

38.965 

Std  Dev 

0.000892 

0.000145 

0.880 

0.093 

0.693 

0IB+M0RB3 

Mean 

0.703752 

0.512869 

18.939 

15.571 

38.799 

Std  Dev 

0.000936 

0.000170 

0.870 

0.093 

0.748 

^Isotopic  variance  is  the  square  of  the  standard  deviation. 
^Mean  and  standard  deviation  based  on  477  samples. 
^Mean  and  standard  deviation  based  on  567  samples. 
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Table  1.3.  Sample  locations  for  the  MORB  data  in  the  OIB+MORB  data  set. 
with  the  number  of  samples  [in  braces]  and  references  indicated. 


■ 

Location 

References 

Atlantic  Ocean  [22] 

2,5 

Pacific  Ocean 

East  Pacific  Rise  [6] 

5,7 

Galapagos  Ridge  [13] 

5,7 

Gorda  Ridge  [8] 

7 

Juan  de  Fuca  Ridge  [6] 

7 

Indian  Ocean  [10] 

1,5 

E  Indian  Ridge  [7] 

4 

SE  Indian  Ridge  [12] 

6 

SW  Indian  Ridge  [6] 

3 

Reference  guide;  [1]  Cohen  and  O'Nions,  1982b:  [21  Cohen  et  ai.  1980;  (31 
Hamelin  and  Allegre,  1985;  [4]  Hamelin  et  al.,  1986;  [5]  Ito  et  ai,  1987;  |6j  Klein 
etal.,  1988;  [7]  White  et  al,  1987. 
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CHAPTER  2 

MATHEMATICAL  AND  STATISTICAL  METHODS  OF  DATA  ANALYSIS 

INTRODUCTION 

When  dealing  with  a  multidimentional  data  set  with  dimension  greater 
than  three,  it  is  impossible  to  visualize  the  shape  of  the  data  in  that  space.  This 
makes  it  difficult  to  choose  "end-members"  for  the  data,  where  end-members  are 
interpreted  as  the  vertices  of  the  smallest  simplex,  with  linear  or  nonlinear  edges, 
that  completely  encloses  all  the  data  points.  Previous  work  using  two- 
dimensional  plots  to  estimate  the  groups  or  end-members  (Zindler  ct  ai.  1982; 
White,  1985;  Zindler  and  Hart,  1986)  can  be  misleading  since  those  plots  are 
projections  of  a  higher-dimensional  shape.  For  this  study,  it  is  possible  to  reduce 
the  dimensionality  of  the  OIB-(-MORB  data  set,  via  principal  component  analysis, 
and  still  retain  its  general  shape,  making  it  possible  to  choose  end-members  in 
three-dimensions. 

For  the  OIB  data  set,  the  data  locations  (oceanic  islands]  are  not  distributed 
evenly  about  the  globe.  This  prompts  the  question  as  to  whether  there  is  any 
relationship  between  location  and  isotopic  signature.  To  address  this,  a  spatial 
correlation  test  (.Mantel,  1967)  is  used  to  test  for  a  correlation  between  the 
geographic  distance  and  the  "isotopic  distance"  between  samples.  In  addition,  a 
count  is  kept  of  the  number  of  times  a  sample's  isotopic  "nearest-neighbor" 
occurs  within  the  same  island  and  within  the  same  geographic  feature. 

Finally,  the  globe  has  been  divided  by  Hart  (1984,  1988)  into  the  islands 
lying  inside  the  DUPAL  belt,  from  2°  S  to  60°  S,  and  those  lying  outside.  To  see 
if  there  is  statistical  justification  for  separating  the  data  into  these  two  different 
populations,  isotopic  nearest-neighbor  discriminant  analysis  is  performed  on  the 
data  set  to  obtain  a  misclassification  error  rate.  The  significance  of  this  error 
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rate  is  based  upon  a  randomization  test  of  Solow  (1990).  While  giving 
promising  results,  the  randomization  test  for  significance  is  inconclusive  because 
spatial  correlation  within  geographic  features  has  not  been  accounted  for.  As  an 
alternative,  discrimination  between  isotopic  signatures  on  the  scale  of  geographic 
features  inside  and  outside  the  DUPAL  belt  is  addressed  graphically. 


PRINCIPAL  COMPONENT  ANALYSIS 

Theory 

Principal  component  analysis  can  be  viewed  as  a  coordinate  system 
transformation,  but  one  that  has  particular  properties.  It  generates  a  new  set  of 
variables,  the  principal  components,  that  are  linear  combinations  of  the  original 
variables: 

5 

=  '  =  ^ . 

7=1 

where  the  Z-s  are  the  principal  components,  the  e^Js  are  the  transformation 
coefficients,  and  the  Xjs  are  the  original  isotope  measurements  (Xj  =  ^^Sr/^^\Sr, 
X2  =  X3  =  206pb/204p5^  X4  =  207pb/204pb^  X5  =  208pt^/2()4pf,) 

The  principal  components  have  the  following  properties: 

(1)  Z^-  and  Zj  are  uncorrelated,  for  all  i,j 

(2)  Variance(Z|)  >  Variance(Z2)  >  .  .  .  >  VarianceCZ^) 

5 

(3)  for  all  /,  X  ^jj  =  1 

7=1 

The  transformation  coefficients  are  the  elements  of  the  unit  eigenvectors  of  the  5 
X  .5  data  covariance  matrix.  Because  the  isotopic  ratios  are  on  different  scales, 
the  data  set  must  be  normalized  in  order  for  all  of  the  isotopes  to  be  treated 
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equally  in  the  analysis.  One  way  to  do  this  is  to  take  each  sample  and  for  every 
isotope  subtract  the  mean  and  divide  by  the  standard  deviation  (Tabic  1 .2): 


where  Xij  is  the  jth  isotopic  ratio  for  the  /th  .  'mple,  etc.  This  method  weights 
the  information  provided  by  all  five  isotopes  equally.  Alternatively,  Allegre  ci 
al.  (1987)  develop  their  owu  empirical  norm,  the  "geologic  norm”,  that  takes 
analytical  errors  into  account  and  is  designed  to  give  equal  weight  to  all  isotopes 
except  2^^Pb/^^‘^Pb,  w  hich  has  the  largest  analytical  error. 


Application  to  the  OIB+MORB  Data  Set 

Because  DMM  (depleted  MORB  mantle]  is  one  of  the  proposed  mantle 
end-member  components,  I  have  chosen  to  do  principal  component  analysis  using 
all  of  the  oceanic  island  data  [477  samples]  plus  a  wide  selection  of  MORB  data 
(90  samples].  The  covariance  matrix  for  the  OIB+MORB  data  .set  and  its 
eigenvectors  and  eigenvalues  are  shown  in  Table  2.1.  The  sum  of  the 
eigenvalues  is  the  trace  of  the  covariance  matrix,  ie.  the  sum  of  the  diagonal 
elements.  This  is  equal  to  5  because  the  diagonal  elements  of  the  covariance 
matrix,  the  scaled  isotope  variances,  are  all  1.  To  find  out  how  much  of  the 
variance  of  the  scaled  data  set  is  accounted  for  by  each  eigenvector,  and  thus  each 
principal  component,  divide  the  corresponding  eigenvalue  by  5.  The  first  three 
principal  components  account  for  97.5%  of  the  variance  of  the  data  set. 

Therefore  it  is  reasonable  to  use  the  three-dimension  principal  component  data 
set  to  select  end-member  components.  This  has  important  implications  for  the 
OIB+MORB  data  set.  In  ^-dimensional  space,  the  polygon  containing  the  fewest 
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vertices  (^t+l  ]  is  a  simplex.  Thus,  the  OIB+MORB  data  set  would  require  six 
end-member  components  to  completely  define  it,  if  it  spanned  the  entire  five¬ 
dimensional  space.  The  fact  that  it  can  be  adequately  represented  in  three- 
dimensions  implies  that  the  OIB-t-MORB  data  set  requires  only  four  end-member 
components. 

A  comparison  of  eigenvalues  and  corresponding  percentages  of  variance 
from  this  study  and  from  Allegre  ct  al.  (1987)  for  OIB-t  MORB  and  OIB  data 
sets  is  presented  in  Table  2.2.  It  should  be  noted  that  the  OIB  eigenvalues  from 
this  study  are  found  using  a  .separate  covariance  matrix  derived  from  the  477 
OIB  samples  alone,  as  is  done  by  Allegre  ct  al.  (1987).  Their  analysis  yielded 
similar  results  for  a  three-dimensional  fit  to  the  data  |01B-t-M0RB;  99.29^  versus 
97. .SC);  OIB:  98.8^7^  versus  97.39^-1.  Part  of  the  small  difference  that  does  exist 
may  be  due  to  the  fact  that  they  used  a  smallei  data  .set  [OIB-i-MORB;  91  samples 
versus  567  samples;  OIB;  53  samples  versus  477],  in  addition  to  the  different 
methods  used  to  scale  the  data. 

I'he  procedure  outlined  above  for  computing  principal  components  is 
compacted  into  matrix  form,  Z  =  F.Y,  with  exact  solutions: 


;z,, ... 

.  .  .  .Z,^ : 

1 

Z''i  . . . 

. . .  .Zzv 
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where  N  =  the  number  of  samples  1567),  the  Y^y's  are  the  normalized  isotopic 
values  and  the  eigenvectors  are  the  rows  of  the  matrix  E.  Three  two- 
dimensional  plots  of  the  first  three  principal  components,  with  general  end- 
member  regions  indicated,  (Figs.  2.1, 2.2  and  2.3)  are  presented  for  comparison 
with  those  of  Allegre  et  al.  (1987)  (Fig.  2.4).  Plotting  the  principal  component 
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values  for  the  samples  versus  each  other  is  the  same  as  plotting  the  projection  of 
the  OIB+MORB  population  onto  its  eigenvector  planes  as  they  have  done.  The 
two  sets  of  plots  are  very  similar,  but  mirror  images  of  each  other.  This  is 
simply  because  the  eigenvectors  used  were  of  opposite  sign,  in  no  way  affecting 
the  validity  of  either  set  of  plots. 

Mantle  End-Member  Components 

In  three-dimensional  space  the  principal  component  data  form  a 
tetrahedron  (Fig.  2.5).  It  should  be  noted  that  the  tetrahedron  is  not  aligned  with 
the  principal  component  axes,  so  two-dimensional  plots  of  the  principal 
component  data  do  not  give  an  exact  indication  of  the  location  of  the  extreme 
points.  End-member  component  values  are  chosen  by  eye  at  the  extretnes  of  tlic 
tetrahedron  using  a  rotating  three-dimensional  plotting  program. 

First,  the  "nonlinear"  end-member  points  are  chosen,  those  that  just  form 
the  vertices  of  the  tetrahedron  (Table  2.3).  These  end-members  are  referred  to 
as  "nonlinear"  because  they  define  the  vertices  of  the  smallest  simplex  enclosing 
the  data  points  which  has  linear  and  nonlinear  edges.  In  geometry,  a  simplex  is 
defined  as  a  polygon  with  planar  faces,  but  I  a?n  extending  this  definition  to 
encompass  polygons  containing  nonplanar  faces  as  well.  The  purpose  of 
choosing  particular  end-member  points  is  to  be  able  to  express  all  of  the  sample 
points  as  a  combination  of  the  four  end-member  components,  for  later  use  in 
spherical  harmonic  expansions.  Though  linear  mixing  is  believed  to  exist 
between  HIMU  and  EMI  (Hart  ct  a/.,  1986)  and  HIMU  and  DMM  (Hart,  1988), 
more  complicated  mixing  arrays  are  probable  amongst  the  other  components. 
Since  no  models  exist  for  the  nonlinear  mixing  arrays,  it  is  easiest  to  represent 
the  sample  points  as  a  linear  combination  of  the  end-member  points.  Thus,  it  is 
necessary  to  find  the  vertices  of  the  smallest  simplex  with  planar  faces  that 
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encloses  as  many  data  points  as  possible;  these  vertices  are  the  "linear "  end- 
members.  These  end-members  are  chosen  by  rotating  the  figure  to  look  at  the 
four  sides  of  the  tetrahedron  edge  on  and  moving  out  the  "nonlinear"  end- 
members  until  the  planar-sided  tetrahedron  defined  by  linear  mixing  expands  to 
contain  as  many  sample  points  as  possible,  without  becoming  overly  extreme 
(Table  2.4).  This  is  an  admittedly  subjective  process,  but  more  accurate  than 
choosing  end-members  using  two  dimensional  plots.  Figures  2.6  -  2.9  show  the 
four  views  normal  to  each  of  the  tetrahedron  faces. 

When  assuming  linear  mixing,  the  simplex  defined  by  the  final  chosen 
"linear"  end-member  points  excludes  only  13  OIB  data  points,  out  of  477,  and  3 
MORE  data  points  out  of  90  (Table  2.5).  compared  to  the  85  OIB  and  49  MORE 
data  points  excluded  when  using  the  "nonlinear"  end-member  values.  The 
excluded  points  will  have  negative  amounts  of  some  of  the  end-members  and  will 
not  be  used  in  spherical  harmonic  expansions. 

The  end-member  values  selected  in  principal  component  space  are 
converted  back  into  normalized  isotope  values  (Tables  2.3  and  2.4)  by 
substituting  zeros  (the  mean  value  for  each  principal  component]  for  the  fourth 
and  fifth  principal  components  in  the  Z  matrix: 


Cii... 

C21 . . . 

.  .  .  .Cl4 

.  .  .  .C24 

^11  ^12^13^.4^15 

-1 

Zii... 
Z21 . . . 

.  .  .  .Zl4 

.  .  .  .Z24 

C31 . . . 

= 

• 

Z3I . . . 

.  . .  .Zv, 

C41  . . . 

. . .  .C44 

-  <^5 1  <^52  ^53  ^54  <^55  - 

0. . .. 

....() 

Xsi  •  ■  • 

- C54_ 

L  0.... 

....()  J 

where  Cj^  is  the  normalized  ^^Sr/^^Sr  ratio  for  the  /th  end-member  component. 

and  so  forth.  There  is  some  error  involved  in  this  process,  but  because  the 
variances  of  the  fourth  and  fifth  principal  components  are  small,  the  error  is 
small.  To  compute  these  errors,  the  entire  OIB+MORB  data  set  is  transfonned 
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into  principal  components;  the  fourth  and  fifth  principal  components  are 
dropped;  the  approximate  normalized  isotope  values  are  computed  as  above;  and 
these  values  are  then  unnormalized  and  compared  to  the  actual  isotope  values. 
The  average  absolute  errors  for  this  transformation  are  fairly  small  compared  to 
the  isotope  standard  deviations  (Table  2.6).  Compared  to  the  range  of  analytical 
errors,  all  of  the  transformation  errors  are  reasonable  except  the  one  for 
206pb/204pb^  which  is  approximately  6x  larger  than  its  analytical  error  (Table 
2.6). 

Finally,  the  samples  are  computed  as  percentages  of  the  four  "linear"  end- 
members; 


Cii  . . . 

. . .  .ChI 

Yi/ 

C21 . . . 

. . .  .C24 

{Plf 

Y2; 

C31 . . . 

. . .  .C34 

'  V 

P2j 

Y3, 

C4I  . . . 

. . .  .C44 

Pis 

Y4; 

C5I . . . 

L  1 .... 

. . .  .C54 

. 1 . 

[P^J 

YS; 

L  1  J 

where  is  the  percentage  of  the  ith  end-member  component  for  the  /th  sample 
and  Y^y  is  the  ith  normalized  isotope  value  for  the  jth  sample.  The  C  matrix  is 

the  normalized  end-member  isotope  value  matrix  computed  from  above  with  an 
additional  row  of  ones.  This  row  of  ones  and  the  one  included  in  the  Y  vector 
define  a  constraint  that  the  sum  of  the  percentages  add  up  to  I .  I’his  is  necessary 
to  provide  useful  positive  results  between  0  and  1  since  the  tetrahedron  is  not  a 
four-component  composition  diagram,  but  resides  in  Euclidean  space.  QR 
decomposition  is  used  to  solve  this  over-determined  system  of  ecjuaiions.  It 
decomposes  the  C  matrix  into  two  matrices;  Q  [orthogonal]  and  R  [upper 
triangular]:  QRp  =  Y,  with  solutions:  =  R'^Q^Y. 
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SPATIAL  CORRELATION  TIESTING 

Methodology 

In  order  to  check  for  spatial  correlation,  a  paired  distance  approach  is 
employed,  as  outlined  in  Mantel  (1967),  using  geographic  and  isotopic  distances. 
The  geographic  distance  used  is  that  of  an  arc  on  a  sphere  connecting  any  two 

sample  locations,  ie.  a  great  circle  distance  (Turcotte  and  Schubert,  1 982).  The 
angle  A^j  between  the  two  locations  /  and  J  on  the  sphere  (Fig.  2.10)  is  given  by; 

Ay  =  COS' '[cos  Qj  cos  0,  +  sin  0y  sin  0,  cos  (tPy-tp/)] 

where  0^  and  (p^-  are  the  colatitude  and  longitude  of  location  /  and  0y  and  tpy  are 
the  colatitude  and  longitude  of  location ./.  The  surface  distance  .v  between  /  and  ./ 
is; 

Sy  =  RAij 

where  R  is  the  radius  of  the  earth  [/?  =  6378.139  km].  The  isotopic  distance  used 
is  the  generalized  Euclidean  distance  in  multidimensions  scaled  by  the  variances 
of  the  isotopic  ratios.  Scaling  by  the  variances  of  the  isotopic  ratios  is  necessary 

to  keep  the  distance  measurement  from  being  dominated  by  the  isotopic  ratio 
with  the  largest  variance,  ^^^Pb/^^'^Pb  (Table  1 .2).  For  any  two  samples  Xj  and 

X^,  the  isotopic  distance  between  them,  d,  is; 


d,y=V(X,-X/V-'(X,-X^) 


where 
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Xi  = 


Xw 

X2/ 

X3/ 

X4/ 

X5J 


is  the  isotope  vector  for  /th  sample  is  the  ^^Sr/^^Sr  ratio  of  the  /th  sample, 
etc.]  and  V  is  the  diagonal  variance  matrix.  A  similar  distance  measurement, 
called  Mahalanobis  distance  (Manly,  1986)  was  considered,  but  not  used  because 
it  utilizes  the  covariance  matrix.  Covariance  is  a  meaningful  measurement  when 
the  data  is  normally  distributed  (elliptical)  in  space.  From  the  three-dimensional 
principal  component  plots  (Figs.  2.5-2.9),  it  is  apparent  that  the  data  set  is  not 
elliptical,  so  covariance  is  a  meaningless  measurement  concerning  the  nature  of 
the  data. 

Next,  the  correlation  between  the  two  distances  for  all  the  samples  is 
calculated.  The  key  to  Mantel's  (1967)  technique  is  to  determine  the  significance 
of  this  observed  correlation  by  creating  random  pairings  of  the  sample  locations 
and  isotopic  signatures,  calculating  the  appropriate  distances,  and  computing 
their  correlation,  thus  constructing  a  distribution  against  which  the  observed 
value  can  be  judged.  This  distribution  is  that  of  the  correlation  under  the  null 
hypothesis  that  the  geographic  distances  are  matched  to  the  isotopic  distances  at 
random. 

Zindler  and  Hart  (1986)  noted  a  relationship  between  the  scale  length  of  a 
geographic  feature  and  the  isotopic  range  of  that  feature.  Basically,  they 
concluded  that  the  largest  isotopic  ranges  occur  in  the  largest  geographic 
features,  while  small  isotopic  ranges  may  occur  in  small  or  large  features.  This 
implies  a  correlation  between  the  within-feature  geographic  distance  and  the 
within-feature  isotopic  distance.  The  paired  distance  correlation  method  outlined 
above  computes  the  correlation  between  geographic  and  isotopic  distances  both 
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within  features  and  between  features.  In  using  this  method,  it  is  possible  that  any 
correlation  within  the  features  may  be  masked  by  a  lack  of  correlation  between 
the  features.  As  an  additonal  check  for  within-feature  correlation,  a  count  is 
kept  of  the  number  of  times  a  sample’s  isotopic  nearest-neighbor  [the  sample  that 
is  the  smallest  isotopic  distance  from  the  sample  in  question]  occurs  within  the 
same  island  and  within  the  same  island  group  [or  island,  if  an  island  is  not  part  of 
a  larger  group].  The  counts  are  performed  both  for  the  observed  data  and  for 
the  random  permutations.  Those  from  the  random  permutations  can  be  used,  as 
before,  to  judge  the  significance  of  the  observed  counts.  The  larger  scale 
geographic  divisions  of  the  data  set  into  island  groups  and  the  remaining  solitary 
islands  (Table  2.7)  will  be  referred  to  from  this  point  on  as  features. 

Application  to  the  OIB  Data  Set 

For  this  application,  the  OIB  data  is  used  since  only  oceanic  island 
interrelationships  are  of  interest.  Two  477  x  477  distance  matrices  are  calculated 
for  the  geographic  and  isotopic  distances  between  samples.  For  the  observed 
data,  the  correlation  between  the  distance  matrices  is  0.1756  and  the  within  island 
and  feature  nearest-neighbor  occurrence  rates  are  61.4%  and  76.7%, 
respectively  (Table  2.8).  The  occurrence  rates  within  islands  and  features  appear 
significant  and  are  confirmed  so  by  randomization,  as  none  of  the  generated 
occurrence  rates  are  as  large  as  the  observed  rates  for  100  permutations  (Table 
2.8).  The  correlation,  on  the  other  hand,  is  small,  but  attains  significance 
compared  to  the  randomization  values  which  are  all  less  than  the  observed  value 
(Table  2.8)  Thus,  both  methods  indicate  that  there  is  spatial  correlation  between 
sample  location  and  isotopic  signature  and  the  correlation  that  exists  between 
samples  within  the  same  geographic  feature  seems  to  dominate. 
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Treating  the  samples  inside  and  outside  the  DUPAL  belt  separately  and 
then  testing  for  spatial  correlation  yields  results  similar  to  those  obtained  with 
the  whole  data  set  (Table  2.8). 

It  is  not  clear  if  all  of  the  spatial  correlation  is  due  to  the  correlation 
within  the  features.  There  may  be  some  additional  spatial  correlation  between 
features.  To  check  this,  the  appropriate  samples  are  averaged  to  get  an  average 
isotopic  signature  and  location  for  each  feature  (Table  2.7).  Using  all  of  the 
features  both  inside  and  outside  the  DUPAL  belt,  the  observed  correlation  is 
0.1584  with  a  significance  level  of  0.13  (there  are  13  permutations,  out  of  100, 
that  have  correlations  higher  than  the  observed  correlation]  (Table  2.8).  Thus,  it 
appears  that  there  is  spatial  correlation  between  features.  However,  if  there  is  a 
distinction  between  features  inside  and  outside  the  DUPAL  belt,  this  distinction 
may  manifest  itself  as  spatial  correlation  when  testing  all  of  the  features  at  once. 
Testing  the  features  inside  and  outside  the  DUPAL  belt  separately  results  in 
correlations  of  0.0685  and  0.2645  with  significance  levels  of  0.95  and  0.51, 
respectively  (Table  2.8).  These  results  indicate  that  there  is  no  significant  spatial 
correlation  between  the  features,  but  that  there  is  a  distinction  between  features 
inside  and  outside  the  DUPAL  belt. 

DISCRIMINANT  ANALYSIS 

Isotopic  Nearest-Neighbor  Discriininant  Analysis 

Methodology.  Without  taking  account  of  spatial  correlation,  the  validity 
of  the  division  of  the  OIB  data  into  samples  inside  and  outside  the  DUPAL  belt  is 
addressed  using  isotopic  nearest-neighbor  as  a  discrimination  rule.  Using  the 
isotopic  distance  measure  outlined  earlier,  a  given  sample's  isotopic  nearest- 
neighbor  is  the  sample  that  is  the  smalle.st  isotopic  distance  away. 
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For  the  discriminant  analysis,  the  assumption  is  made  that  the  selected 
sample's  location  is  unknown,  so  it  is  assigned  the  location  of  its  isotopic  nearest- 
neighbor.  This  assigned  location  is  compared  to  the  actual  location;  if  they  are 
different,  it  is  a  misclassification.  A  count  is  kept  of  the  number  of 
misclassifications  to  calculate  an  error  rate. 

Solow  (1990)  proposes  a  randomization  technique  forjudging  the 
estimated  misclassification  probability  or  error  rate.  The  importance  of  the 
misclassification  error  rate  is  to  test  the  null  hypothesis  that  there  is  no 
difference  between  the  samples  inside  and  outside  the  DUPAL  belt.  This  is  a 
trivial  matter  if  the  sampling  distribution  of  the  error  rate  under  the  null 
hypothesis  is  known,  but  in  this  case  it  is  not.  A  simple  but  effective  way  to 
judge  the  significance  of  the  observed  error  rate  is  to  construct  a  randomization 
distribution  under  the  null  hypothesis  that  the  pairing  of  i.sotopic  signatures  and 
locations  inside  or  outside  the  DUPAL  bell  occurs  by  chance.  Applying  the 
.  randomization  technique  to  the  data,  the  samples  retain  their  isotopic  signature, 
so  their  isotopic  nearest-neighbor  remains  the  same,  but  they  are  randomly 
assigned  to  locations  inside  and  outside  the  belt.  The  discriminant  analysis  is 
done,  as  described  above,  with  this  new  randomly  constructed  data  set  to  get  its 
misclassification  error  rate.  Then  the  process  is  repeated  to  construct  the 
distribution. 

Application  to  the  OIB  data  set.  For  the  OIB  data  set,  the  ob.served 
misclassification  rate  is  7.3%  and  the  randomization  error  rate  ranges  from 
35.2%  to  53.7%.  Superficially,  it  appears  that  describing  the  data  as  two 
populations  residing  inside  and  outside  the  DUPAL  belt  is  viable.  However,  the 
within-feature  spatial  correlation  has  not  been  accounted  for  in  this  analysis.  If 
76.7%  of  the  time,  a  sample's  isotopic  nearest-neighbor  is  located  within  the 
same  geographic  feature,  then  it  seems  obvious  that  the  misclassification  error 
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rate  would  be  small.  The  observed  error  rate  itself  is  not  incorrect,  but  the 
randomization  distribution  of  error  rates  against  which  it  is  being  judged  is 
incorrect.  In  order  for  the  significance  of  the  observed  error  rate  to  be  properly 
judged,  the  spatial  correlation  must  be  preserved  in  the  randomization  process. 

In  this  case,  preserving  the  spatial  correlation  is  too  complicated  to  pursue  when 
other  methods  may  provide  the  desired  information. 

Graphical  Discrimination  of  Geographic  Features 

As  shown  earlier,  the  correlation  between  isotopic  distance  and  the 
geographic  distance  within  features  is  very  strong.  A  way  around  this  spatial 
correlation  is  to  look  for  differences  between  populations  inside  and  outside  the 
DUPAL  belt  on  the  feature  level.  The  averaged  isotopic  values  for  the  features 
(Table  2.7)  are  scaled  by  the  mean  ana  standard  deviation  of  the  isotopes  derived 
from  the  entire  OIB+MORB  data  set  (Table  1.2)  and  expressed  in  terms  of 
principal  components  using  the  eigenvectors  of  the  OIB+MORB  correlation 
matrix  (Table  2.1). 

The  first  three  principal  components  are  plotted  to  look  for  differences  in 
features  inside  and  outside  the  DUPAL  belt,  with  the  general  direction  of  the 
end-members  indicated  (Figs.  2.1 1-2.13).  On  all  of  the  plots,  but  especially  Z3 
versus  Z2,  most  of  the  features  outside  the  belt  cluster  in  a  band  between  DMM 

and  HIMU,  with  the  exception  of  the  Hawaiian  islands  [the  Koolau  volcanics  on 
Oahu  show  a  strong  EMI  signature  (Hart,  1988)],  Shimada  Seamount  [which  aiso 
has  an  EMI  signature  (Hart,  1988)],  and  the  Azores  [Sao  Miguel  has  a  strong 
EMII  signature  (Hart,  1988)].  Essentially,  the  features  outside  the  DUPAL  belt, 
with  few  exceptions,  occupy  only  part  of  the  available  isotopic  space,  while 
features  inside  the  belt  occupy  all  of  the  available  isotopic  space,  including  some 
overlap  with  features  outside.  This  is  essentially  the  relationship  found  by  Hart 
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(1988),  not  that  the  two  populations  are  totally  separated,  but  that  one  population 
contains  isotopic  signatures  that  the  other  does  not.  It  is  important  that  this  two 
population  distinction  is  still  valid  on  the  feature  level.  Since  it  is  still  apparent 
at  this  larger  scale  [not  just  sample  to  sample]  the  geochemical  signatures  of  the 
oceanic  island  basalts  do  have  a  long  wavelength  component  to  them,  making  it 
feasible  to  attempt  to  quantify  these  signatures  using  spherical  harmonic 
expansions. 

In  addition  to  this  graphical  presentation,  the  discrimination  analysis  can 
also  be  done  on  the  feature  level,  but  the  variances  of  the  isotopes  within  each 
feature  must  be  accounted  for  in  some  way. 

SUMMARY 

Mathematical  and  statistical  methods  to  explore  and  characterize  the  01 B 
and  MORB  data  reveal  these  main  points: 

•  OIB+MORB  data  require  only  four  mantle  end-member  components  to 
completely  span  the  range  of  known  isotopic  values. 

•  Choosing  the  mantle  end-member  components  can  be  made  easier  (and 
more  accurate)  with  the  use  of  principal  component  analysis. 

•  Within  geographic  features,  there  is  a  significant  correlation  between 
location  and  isotopic  signature,  but  between  geographic  features,  there 


IS  not. 
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•  Graphical  discrimination  of  geographic  features  shows  that  the 
distinction  between  islands  inside  and  outside  the  DUPAL  belt  is 
viable. 

•  The  existence  of  the  DUPAL  anomaly  on  the  feature  level  indicates 
that  the  anomaly  has  a  long  wavelength  component  to  it. 
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Table  2.1.  Covariance  matrix^  of  the  five  isotopes  with  its  eigenvectors  and 
eigenvalues. 


Covariance  Matrix 


Isotope 

Yl 

Y2 

Y3 

Y4 

Ys 

Yi 

1.000000 

-0.796442 

-0.273004 

-0.019107 

0.061599 

Y2 

1.000000 

0.054987 

-0.170370 

-0.295078 

Y3 

1.000000 

0.901205 

0.894577 

Y4 

1 .000000 

0.901429 

Y5 

1  .OOOOOO 

Eigenvector  Matrix 

isotope 

I 

II 

III 

IV 

V 

Yl 

-0.017647 

0.699432 

0.682362 

-0.174033 

-0.120738 

Y2 

-0.122196 

-0.683457 

0.679315 

-0.179201 

0.156120 

Y3 

0.565079 

-0.195210 

0.019019 

-0.235936 

-0.765867 

Y4 

0.574974 

-0.006763 

A  A  AAOn 

A  "Tr'-»Ar%n 

t  ,  O 

0.200145 

Ys 

0.578661 

0.074352 

-0.102014 

-0.561051 

0.578307 

Eigenvalues 

2.830 

1.861 

0.183 

0.091 

0.035 

Percentage  of  total  variance  accounted  for  b>  each  eigenvector 

56.6  37.2  3.7  1.8  0.7 


*Only  the  upper  half  of  the  covariance  matrix  is  shown  since  it  is  symmetric. 

All  eigenvector  values  are  rounded  to  six  decimal  places  from  the  fourteen 
decimal  accuracy  used  in  the  calculations. 

Covariance  matrix  is  calculated  using  477  OIB  and  90  MORB  samples. 
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Table  2.2.  Comparison  of  eigenvalues  and  percentages  of  variance  accounted  for 
by  the  corresponding  eigenvectors  from  this  study  and  from  Allegre  ct  al. 

(1987)’  for  OIB+MORB  and  OIB  data  sets. 


OIB+MORB 

I 

II 

III 

IV 

V 

2 

2.830 

1.861 

0.183 

0.091 

0.035 

[56.6%] 

[37.2%] 

[3.7%] 

[1.8%] 

[0.7%,] 

^Allegre 

3.20 

1.61 

0.15 

0.03 

0.01 

et  al. 

[64.0%] 

[32.2%] 

[3.0%] 

|0.6%] 

[0.2%  1 

OIB 

I 

II 

III 

IV 

V 

4 

3.047 

1 .568 

0.249 

0.099 

0.037 

[60.9%.] 

[31.4%] 

[5.0%] 

(2.0%! 

|().7%| 

'"'Allegre 

2.85 

1.87 

0.22 

0.05 

0.01 

et  al. 

[57.0%] 

[37.4%] 

[4.4%] 

[1.0%] 

|0.2%4 

’Eigenvalues  from  Allegre  et  al.  (1987)  are  converted  to  scaled  eigenvalues  that 
add  up  to  5  for  comparison  with  eigenvalues  from  this  study. 

Percentages  of  variance  accounted  for  by  the  corresponding  eigenvectors  are 
indicated  in  parentheses. 

^Based  on  567  samples. 

^Based  on  91  samples. 

'^Based  on  477  samples. 

"'Based  on  53  samples. 
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Table  2.3.  "Nonlinear"  end-member  component  values  in  principal  component 
space  and  the  transformed  values  in  isotope  space. 


End-Members  in  Principal  Component  Space 


Z1 

Z2 

Z3 

EMI 

-2.0 

3.6 

-1.3 

EMU 

1.0 

4.0 

2.2 

HIMC 

5.0 

-1.3 

-0.25 

DMM 

-3.75 

-2.9 

0.4 

End-Members 

in  Isotope  Space 

XI 

X2 

X3 

X4 

X5 

:mi 

0.705311 

0.512343 

17.322 

15.431 

38.232 

;mii 

0.707759 

0.512638 

18.788 

15.673 

39.287 

IIMU 

0.702659 

0.512887 

21.615 

15.833 

40.91  1 

►MM 

0.702171 

0.513329 

17.594 

15.381 

36.983 
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1'able  2.4.  "Linear"  end-member  component  value.s,  based  upon  linear  nii.xing, 
in  princi[ial  component  space  and  the  transformed  values  in  isotope  space. 


End-Members  in  Principal  Component  Space 


Z1 

Z2 

Z3 

EMI 

-2.4 

3.6 

-1.6 

EM  1 1 

1.8 

4.5 

2.6 

HIMU 

6.0 

-1.9 

-0.6 

DMM 

-4.3 

-3.7 

0.35 

End-Members  in  Isotope  Space 


XI  X2  X3  X4  X5 


EMI 

0.705126 

0.512316 

17.121 

15.403 

38.082 

EMII 

0.708329 

0.512609 

19.103 

15.724 

39.630 

HIMU 

0.702026 

0.512896 

22.203 

15.879 

41.337 

DMM 

0.701624 

0.513428 

17.459 

15.351 

36.704 
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Table  2.5.  Samples  excluded  from  linear  mixing  tetrahedral  volume.* 


Location 

Sample  Number 

Row  Number^ 

Azores,  Sao  Miguel 

SMID 

32 

SM49 

36 

Galapagos 

E35 

173 

Gough 

10 

175 

Hawaii 

69Tan2 

200 

Kerguelen  Plateau 

DR02/12 

279 

DR08 

282 

747c-12r-4-45-46 

292 

747c- 16r-2-8 1-84 

294 

Marquesas 

uapl  1 

329 

Pitcairn,  Pulwana 

642 

370 

St.  Helena 

37 

469 

237 

482 

Atlantic  Ocean 

AD3-3 

535 

SW  Indian  Ridge 

D1 

536 

D5 

539 

*OIB  samples  excluded  from  the  volume  will  not  be  used  in  spherical  harmonic 
expansions. 

^Indicates  row  number  of  the  data  set  included  in  Appendix  A. 


44 


Table  2.6.  Average  absolute  errors  in  transforming  three-dimensional  principal 
component  data  into  five-dimensional  isotope  data,  with  their  ratio  to  isotope 
standard  deviations  and  comparison  to  analytical  errors. 


• 

XI 

X2 

X3 

X4 

X5 

Average  absolute  error  [for  567  samples] 

0.000041  0.000008  0.111  0.016  0.111 

Ratio  of  average  absolute  error  to  standard  deviation^ 

0.043802  0.047124  0.128  0.170  0.149 


Absolute  Error  Percentage  Range^  Analytical  Error  Range 


XI 

0.00580  to  0.00584% 

0.003  to  0.01% 

X2 

0.001558  to  0.001561% 

0.003  to  0.01% 

X3 

0.225  to  0.330%/amu 

0.03  to  0.05%/amu 

X4 

0.0337  to  0.0347%/amu 

0.03  to  0.05%/amu 

X5 

0.068  to  0.075%/amu 

0.03  to  0.05%7amu 

*  Isotopic  standard  deviations  for  the  OIB-(-MORB  data  set  are  indicated  in  Table 
1.2. 

“Absolute  error  percentage  ranges  are  calculated  using  the  average  absolute 
errors  and  the  ranges  of  the  isotopes  in  the  OIB-t-MORB  data  set; 

XI  0.702290  to  0.707400 
X2  0.512376  to  0.513290 

X3  16.943  to  21.755 

X4  15.406  to  15.862 

X5  37.235  to  40.619 


Table  2.7.  Average  isotopic  signatures  and  locations  of  the  geographic  features  (island  croups,  islands,  ridces,  seamounts) 
represented  in  the  OIB  data  set  with  the  number  of  samples  (in  braces). 


Crozet  Islands  [9]  0.703997  0.512849  18.929  15.587  39.037  -46.45  52.00 
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Table  2.7.  Continued. 
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Table  2.8.  Correlations  between  geographic  and  isotopic  distance  matrices  and 
island/feature  isotopic  nearest-neighbor  occurrence  rates  for  all  samples  in  the 
OIB  data  set  and  samples  inside  and  outside  the  DUPAL  l>"h  rorrelations  for  all 
the  geographic  features  and  those  inside  and  outside  the  DUPAL  belt  are  also 
given. 


Correlation 

Randomization  Range^ 

OIB 

0.1756 

-0.0077  to  0.0224 

Inside  DUPAL 

0.0641 

-0.0035  to  0.0343 

Outside  DUPAL 

0.6142 

-0.0045  to  0.0699 

Island  Occurrence  Rate 

Randomization  Range^ 

OIB 

61.4% 

1.0%  to  3.6% 

Inside  DUPAL 

63.6% 

1 .6%  to  6.2% 

Outside  DUPAL 

47.1% 

1.7%  to  11.0% 

Feature  Occurrence  Rate 

Randomization  Range^ 

OIB 

76.7% 

3.4%  to  10.0% 

Inside  DUPAL 

71.5% 

4.2%  to  11.1% 

Outside  DUPAL 

75.6% 

19.2%  to  35.5% 

Correlation  Randomization  Range^  Significance 

Level 


Features 

0.1584 

0.0369  to  0.2075 

0.13 

Inside  DUPAL 

0.0685 

0.0350  to  0.2807 

0.95 

Outside  DUPAL 

0.2645 

0.0649  to  0.6602 

0.51 

1 


Randomization  ranges  based  upon  100  random  permutations. 
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Z1 


Fig.  2.1.  Plot  of  the  second  principal  component  versus  the  first  principal 
component  for  the  OIB+MORB  data  set.  Symbols:  x  =  MORB  data,  open  circle 
=  OIB  samples  inside  the  DUPAL  belt,  black  diamond  =  samples  outside  the 
DUPAL  belt.  General  mantle  end-member  component  regions  are  indicated. 
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Fig.  2.3.  Plot  of  the  third  principal  component  versus  the  second  principal 
component  for  the  OIB+MORB  data  set.  Symbols:  x  =  MORB  data,  open  circle 
=  OIB  samples  inside  the  DUPAL  belt,  black  diamond  =  samples  outside  the 
DUPAL  belt.  General  mantle  end-rr.ember  component  regions  are  indicated. 
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VI 


Fig.  2.4.  Plots  of  the  second  principal  [V2]  component  versus  the  first  fVl]  and 
the  third  principal  component  [V3]  versus  the  first  fVl]  for  a  smaller 
OIB+MORB  data  set  from  an  analysis  done  by  Allegre  et  al.  (1987).  These  plots 
are  the  mirror  images  of  the  ones  done  for  this  analysis  because  the  chosen 
eigenvectors  for  the  two  analyses  are  of  opposite  sign. 


Fig.  2.5.  Three-dimensional  view  of  the  OIB-i-MORB  principal  component  data. 
Axes:  X  =  ZI,  Y  =  Z2,  Z  =  Z3.  Symbols  for  the  end-member  components:  +  = 
EMI,  X  =  EMU,  diamond  =  HIMU,  square  =  DMM. 
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Geographic 

North 


Fig.  2.10.  Geometry  for  determining  the  surface  distance  s  between  locations  I 
and  J  on  the  globe,  where  0  and  (p  are  colatitude  and  longitude  and  A  is  the  angle 
between  the  two  locations  taken  from  the  center  of  the  Earth.  From  Turcotte 
and  Schubert  (1982). 
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Z1 


Fig.  2.1 1.  Plot  of  the  second  principal  component  versus  the  first  principal 
component  for  the  36  geographic  features.  Symbols:  open  circle  =  features 
inside  the  DUPAL  belt,  black  diamond  =  features  outside  the  DUPAL  belt. 
Labeled  points:  1  =  Hawaiian  Islands,  2  =  Shimada  Seamount,  3  =  Azores.  The 
general  directions  of  the  mantle  end-member  component  regions  are  indicated. 
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Fig.  2.1Z  Plot  of  the  third  principal  component  versus  the  first  principal 
component  for  the  36  geographic  features.  Symbols:  open  circle  =  features 
inside  the  DUPAL  belt,  black  diamond  =  features  outside  the  DUPAL  belt. 
Labeled  points:  1  =  Hawaiian  Islands,  2  =  Shimada  Seamount,  3  =  Azores.  The 
general  directions  of  the  mantle  end-member  component  regions  are  indicated. 
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Fig.  2.13.  Plot  of  the  third  principal  component  versus  the  second  principal 
component  for  the  36  geographic  features.  Symbols:  open  circle  =  features 
inside  the  DUPAL  belt,  black  diamond  =  features  outside  the  DUPAL  belt. 
Labeled  points:  1  =  Hawaiian  Islands,  2  =  Shimada  Seamount,  3  =  Azores.  The 
general  directions  of  the  mantle  end-member  component  regions  are  indicated. 
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CHAPTER  3 

SPHERICAL  HARMONIC  REPRESENTATION  OF  ISOTOPIC 
SIGNATURES:  THE  CONTINUOUS  LAYER  MODEL 

INTRODUCTION 

Hart  (1984)  contoured  world  maps  of  OIB  isotope  data  for  his  three 
DUPAL  anomaly  criteria  (ASr  >  40;  A7/4  >  3;  A8/4  >  40],  These  maps  show  a 
concentrated  band  spanning  approximately  60°  of  latitude,  centered  on  30°-40°S, 
with  pronounced  highs  for  the  anomaly  criteria  in  a  region  from  the  South 
Atlantic  to  the  Indian  Ocean  [ASr,  A7/4,  A8/4]  and  in  the  central  Pacific  [ASr, 
A8/4].  Qualitatively,  Hart  (1984,  1988)  believes  this  geochemical  anomaly 
correlates  with  other  geophysical  anomalies:  the  slab-corrected  geoid  (Hager, 
1984),  deep  mantle  P-wave  tomography  maps  (Dziewonski,  1984),  slow  P-wave 
regions  at  the  core/mantle  boundary  (Creager  and  Jordan,  1986)  and  equatorial 
anomalies  in  the  core  (Le  Mouel  et  ai,  1985).  These  geophysical  anomaly 
patterns  are  typically  expanded  in  terms  of  spherical  harmonics,  therefore  any 
attempt  to  make  a  quantitative  comparison  between  geochemical  and  geophysical 
patterns  requires  expanding  the  geochemistry  data  in  spherical  harmonics  as 
well. 

Expansion  of  the  geochemistry  data  is  approached  in  two  ways,  based 
upon  an  assumed  geometry  for  the  OIB  geochemical  reservoir.  The  first 
approach,  the  "continuous  layer  model"  discussed  in  this  chapter,  assumes  that 
the  OIB  reservoir  is  a  continuous  layer  [not  ruling  out  heterogeneities  within  this 
layer]  and  tries  to  reconstruct  this  layer.  Plumes  from  this  layer  only  sample  the 
continuous  geochemical  "function"  in  discrete  locations.  With  the  geochemistry 
"function"  unknown,  the  spherical  harmonic  coefficients  must  be  solved  for 
using  least  squares,  singular  value  decomposition  or  a  similar  method  that  will 
approximate  the  values  of  the  geochemistry  "function"  where  there  is  no  data. 
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'I'hc  second  approach,  the  "delta-function  model"  discussed  in  Chapter  4,  assumes 
that  the  OIB  reservoir  is  composed  of  a  series  of  point  sources,  each  feeding  a 
separate  plume.  In  this  case,  the  geochemistry  "function"  is  known  and  can  be 
represented  as  a  series  of  delta-functions.  The  spherical  harmonic  coefficients 
can  be  solved  for  directly  with  the  simplification  from  integration  to  summation 
allowed  by  the  delta-function  approximation. 

The  continuous  layer  model  and  the  delta-function  model  are  not  meant  to 
suggest  two  end-member  possibilities  for  OIB  source  geometry.  Rather,  the 
delta-function  model  can  be  regarded  as  an  approximation  of  the  continuous 
layer  model  that  gives  a  mathematically  robust  solution  for  the  spherical 
harmonic  coefficients.  In  regard  to  the  oceanic  crust  model  of  Hofmann  and 
White  (1982),  the  continuous  layer  model  corresponds  to  the  accumulated  layer 
of  subducted  oceanic  crust,  with  the  plume-forming  instabilities  occurring  at 
discrete  locations  within  this  layer.  The  delta-function  model  can  also  be 
reconciled  with  the  accumulated  layer  model,  with  the  stipulation  that  discrete 
pockets  [point  sources]  within  this  layer  form  and  feed  individual  instabilities. 

For  the  purposes  of  minimizing  small  scale  variations  [ie.  variations 
within  a  single  island  or  island  chain]  in  the  geochemistry  "function”  that  cannot 
be  accurately  represented  with  the  incomplete  elobal  data  coverage,  this  spherical 
harmonic  study  is  based  on  the  averaged  isotopic  signatures  of  the  36  geographic 
features  (Table  2.6).  These  average  isotopic  values  are  converted  to  mantle-end 
member  component  percentages  (Table  3.1),  as  outlined  in  Chapter  2,  to  form 
the  data  matrices  used  in  the  expansions. 


SPHERICAL  HARMONIC  BASICS 

Spherical  harmonics.  Y'/'(0,(p),  are  a  set  of  orthonormal  functions  over  the 


unit  sphere: 


65 


Yr(e,<p)  =  17(^056)  e.'”<P 

V  47i(/+w)! 

where  /  is  the  degree  of  the  expansion,  m  is  the  order  of  the  expansion,  0  is 
colatitude  [0  =  71/2  -  latitude;  O<0<7r]  and  (p  is  longitude  (-7C<(p<7i).  The 
functions  form  a  complete  set  of  orthogonal  functions  in  the  index  m  on  the 
interval  -7r<(p<7r  and  the  associated  Legendre  polynomials  P^(cos0)  form  a 
similar  set  in  the  index  /  for  each  m  value  on  the  interval  -l^os0<l  (Jackson, 
1975).  Therefore  their  product  forms  a  complete  orthogonal  set  on  the  surface 
of  the  unit  sphere  in  the  two  indices  l,m.  The  spherical  harmonic  functions  used 
in  this  analysis  are  normalized  by  the  square  root  term  so  that  their  integrated 
square  over  the  sphere  is  unity  fin  most  geophysics  applications,  the  functions  are 
normalized  so  that  the  integrated  square  over  the  sphere  is  471]: 

I  dtp  I  r/(cos0)\7'(0,(p)‘Yf(0,(p)  =  6/76,„v. 


where  the  asterisk  denotes  complex  conjugation. 

Any  function  /(0,cp)  can  be  expanded  in  spherical  harmonics: 

/■(e.-p)  =  i  i 

/=()  m=-t 


where  L  is  the  maximum  degree  of  the  expansion  and  CP  are  complex  spherical 
harmonic  coefficients.  Written  in  a  more  explicit  form,  the  equation  becomes: 


= t  i 


M)  tn=^() 


/(2/+l)(7-m7 

V  47r(/+/;0! 


P^(cos0)  [  cos/?/(p  +  By^sin»j(pj 
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where  A/"  and  B"'  are  real  spherical  harmonic  coefficients.  When  expanding  a 
function  from  degrees  0  to  L,  the  number  of  coefficients  that  need  to  be 

12/+! 

calculated  is:  There  are  actually  an  additional  [L+1  ]  coefficients 

involved,  but  for  m  -  0,  sinmtp  =  0,  so  =  0.  It  is  important  to  realize  that 
only  having  36  features  limits  the  possible  spherical  harmonic  expansion  to 
degree  5,  in  order  to  avoid  a  purely  underdetermined  problem. 


MANTLE  END-MEMBER  COMPONEN  I'S 
When  attempting  to  use  inverse  methods  to  solve  for  the  htirmonic 
coefficients  of  an  unknown  function,  careful  attention  must  be  paid  to  the 
variation  of  the  data  as  a  function  of  distance  to  avoid  the  problem  of  aliasing. 
For  a  simple  two-dimensional  case,  aliasing  occurs  if  the  sampling  interval  is 
longer  than  half  the  shortest  wavelength  of  the  function  sampled,  causing  the 
sampled  points  to  show  a  periodicity  that  docs  not  exist  in  the  original  data.  The 
minimum  distance  between  any  two  geographic  features  in  the  OIB  feature  data 
set  is  33.396  km,  but  the  distance  between  features  is  not  constant.  Plots  of  data 
variation  versus  distance  between  data  locations  make  it  possible  to  select  a 
mininum  sampling  distance  based  on  the  shortest  distance  required  to  get  the 
maximum  data  variation.  This  minimum  sampling  distance  then  controls  the 
minimum  degree  to  which  the  data  must  be  expanded  in  order  to  adecpiately 
represent  the  data  in  spherical  harmonics  without  aliasing.  The  relationship 
between  wavelength  and  degree  is; 
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where  X  is  the  wavelength  [>i  =  2*(sampling  distance)],  R  is  the  radius  of  the 
earth  [/?  =  6378.139  km]  and  I  is  degree.  Solving  for  degree  in  terms  of 
wavelength: 


2 


Variation-Distance  Relationships 

For  variation-distance  relationships,  the  distance  measure  is  the  angle 

[in  degrees]  from  the  center  of  the  earth  between  any  two  locations  /  and ./  J.see 

Chapter  2]  and  the  variation  measure  is  the  absolute  value  of  the  difference 
between  the  mantle  component  percentages  at  those  locations.  The  angle  A^j  can 

be  transformed  into  a  great  circle  distance  in  km  by  converting  A^j  to  radians  and 

multiplying  by  the  radius  of  the  earth  R. 

Plots  of  absolute  difference  versus  angle  for  the  four  mantle  components 
(Figs.  3.1 -3.4)  show  the  maximum  variation  in  the  components  occurring  on 
very  short  distance  scales  for  the  EMI  and  HIMU  components  and  moderate 
distance  scales  for  the  EMU  and  DMM  components.  Based  upon  these  plots,  the 
minimum  sampling  distances  (in  degrees]  are  ~  14.5°  for  EMI  and  HIMU,  ~  39° 
for  EMU  and  ~  57°  for  DMM.  These  correspond  to  expansions  out  to  degrees 
12,  4  and  3,  respectively.  For  the  current  problem,  the  EMU  and  DMM  data  sets 
can  be  expanded  in  spherical  harmonics  as  they  are,  but  the  EMI  and  HIMU  data 
.sets  require  some  additional  manipulation. 

Variation  Reduction  by  Categorizing  Features 

Separation  of  the  geographic  features  into  populations  located  inside  and 
outside  the  DUPAL  belt  |2°S  to  60°S]  does  not  result  in  two  distinct  isotopic 
populations  [Chapter  2j.  Essentially,  one  population  (outside  the  belt]  defines  a 
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small  field  in  isotopic  space,  while  the  other  population  [inside  the  belt)  defines  a 
larger  field  that  overlaps  with  the  smaller  field  (Fig.  3.5).  A  possible  source  of 
the  large,  small-scale  isotopic  variation  exhibited  by  the  EMI  and  HIMU  data  sets- 
is  the  juxtaposition,  due  to  the  overlap  in  isotope  space,  of  features  having  a 
strong  DUPAL  signature  next  to  those  that  do  not.  If  it  is  possible  to  separate 
DU  PAL-type  features  [those  features  showing  a  strong  DUPAL  signature]  from 
DMM-type  features,  this  separation  might  reduce  the  small-scale  variation  within 
these  two  populations  and  thus  reduce  the  degree  to  which  the  population  data 
must  be  expanded. 

Since  the  goal  is  to  separate  DUPAL-type  features  from  DMM-type 
features,  a  logical  starting  place  is  to  look  at  the  spatial  distribution  of  different 
percentage  categories  of  the  DMM  component  in  three-dimensional  principtil 
component  space  (Fig.  3.6).  Six  DMM  percentage  categories  |<I()%,  10-20%, 
20-30%,  30-40%,  40-50%,  >50%]  can  be  distinguished  as  six  separate  point 
groupings.  Most  striking  is  a  large  spatial  separation  that  occurs  within  the  30- 
40%  category  for  a  small  percentage  difference  [Louisville  -  31.84%,  Balleny  - 
32.17%,  Cocos  -  38.53%].  This  is  a  reasonable  place  to  separate  the  DUPAL- 
type  features  from  the  DMM-type  features,  with  a  boundary  value  of  32% 

DMM,  for  simplicity.  The  resulting  27  DUPAL-type  features  and  9  DMM-type 
features,  with  their  percentage  of  the  DMM  component  are  listed  in  Table  3.2. 

There  are  too  few  DMM-ty{)e  features  to  draw  any  conclusions  from  plots 
of  absolute  difference  versus  angle.  For  the  DUPAL-type  features,  plots  of 
absolute  difference  versus  angle  of  the  DUPAL  components  [EMI,  EMU  and 
HIMU]  show  no  reduction  in  the  small-scale  variation,  while  that  of  DMM  does, 
with  an  increase  in  sampling  distance  from  ~  57°  to  89°  (Figs.  3.7-3.10).  In 
retrospect,  this  is  an  obvious  result  of  the  artificial  separation  performed.  3'hc 
percentage  categories  are  basically  parttllcl  slices  through  the  tetrahedron  th;ii 
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move  from  a  broad  base  of  lower  percentages  to  a  peak  of  high  percentages 
approaching  an  end-member  component  apex  on  the  tetrahedron  [like  a  ternary 
diagram].  It  is  true  that  these  slices  can  separate  DMM-type  features  from 
DUPAL-type  features,  but  only  the  variation  of  the  DMM  percentages  are 
reduced.  To  reduce  the  variation  of  the  individual  DUPAL  components  using 
this  method,  EMI-type  features  would  have  to  be  distinguished  from  non-EMI- 
type  features,  etc.  This  would  generate  four  different,  though  overlapping,  sets 
of  features  to  use  to  characterize  the  four  different  components.  Manipulation  of 
the  data  set  in  this  way  is  not  desirable,  so  another  method  must  be  pursued  in  tlic 
attempt  to  reduce  small-scale  data  variation. 

Variation  Reduction  by  Filtering 

Another  method  to  reduce  small-scale  variation  [and  hopefully  enhance 
any  long  wavelength  component]  is  to  filter  the  data  set  in  some  way.  Here,  a 
simple  circular  filter,  of  fixed  radius,  is  applied  to  each  feature  location.  The 
new  data  values  assigned  to  that  feature  location  are  the  means  of  the  mantle 
component  percentages  of  the  feature  locations  that  fall  within  the  circle.  To 
ensure  that  there  are  always  at  least  two  features  falling  within  the  circle,  the 
radius  of  this  circle  is  determined  by  the  longest  distance  to  the  nearest  feature 
location.  Nunivak  Island  is  the  most  isolated  feature  v'ith  the  nearest  feature 
being  the  Hawaiian  Islands  at  an  angular  distance  [from  the  center  of  the  earth]  of 
40.86°.  The  circle  radius  is  then  40.9°,  for  simplicity. 

Plots  of  absolute  difference  versus  angle  for  the  filtered  data  set  yield 
interesting  results  (Figs.  3.1 1-3.14).  All  of  the  mantle  component  data  sets  shou 
a  reduction  in  small-scale  variation,  except  EMIl,  which  shows  an  increase  in 
\  ariation,  with  a  decrease  in  angular  sampling  distance  from  ~  39°  to  27° 
K'xpansions  to  degrees  4  and  7,  respectively |.  I'he  remaining  plots  show  an 


increase  in  angular  sampling  distance  from  ~  14.5°  to  37°  [expansions  to  degrees 
12  and  5,  respectively]  for  EMI,  an  increase  from  ~  57°  to  102°  (expansions  to 
degrees  3  and  2,  respectively]  for  DMM  and  a  dramatic  increase  from  ~  14.5°  to 
83°  [expansions  to  degrees  12  and  2,  respectively]  for  HIMU.  Now  that  the 
small-scale  variation  has  been  significantly  reduced  by  filtering,  the  filtered  EMI 
and  HIMU  data  sets  can  also  be  expanded  in  spherical  harmonics. 

INSIGHTS  FROM  GEOPHYSICAL  DATA 
It  is  unclear  how  accurate  the  spherical  harmonic  expansions  of  the  01 E 
feature  data  set  will  be  due  to  the  limited  global  coverage  and  the  highly  variable 
nature  of  the  data.  In  an  attempt  to  address  these  problems,  three  geophysical 
data  sets,  with  different  variance  characteristics,  are  constructed  with  the  same 
limited  coverage  to  provide  a  sort  of  control  set  against  which  qualitative 
comparisons  can  be  made.  Geoid,  gravity  and  gravity  gradient  anomalies  are  the 
chosen  geophysical  measures  because  their  coefficients  are  well  known  and  they 
form  a  kind  of  continuum  from  the  long  wavelength  ]low  degree]  dominance  in 
the  geoid  signature  to  the  short  wavelength  (high  degree]  dominance  in  the 
'.ravity  gradient  signature  (Fig.  3.15).  Techniques  applied  to  the  mantle 
)mponent  data,  to  solve  for  the  spherical  harmonic  coefficients,  are  also  applied 
these  constructed  data  sets  to  see  how  closely  the  actual  geophysical 
fficients  can  be  approximated. 

truction  of  (Jeophysical  Data  Sets 
The  gravitational  potential  V^,  in  spherical  harmonics  as  a  function  of 
listance  r,  is  given  by: 
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V/  .  GK  X  Ff  (cos0)[a7’cosw(P  +  BT^sinmcp]! 

R  \'  1=2  '  n.=o  V  4Tc(/+m)!  I 

where  G  is  the  gravitational  constant  [G  =  6.6726x10'^  ^  m^/kg»s^J,  M  is  the 
mass  of  the  earth  [M  =  5.973x10^'^  kg]  and  R  is  the  radius  of  the  earth  in  meters 
(Stacey,  1977).  The  gravitational  potential  anomaly  [<5^  =  ^observed  ' 

^""theoretical^ 

.  QM- 1  y  '  y  pw  [AT’cosmtp  +  Bf'sinmtp] ! 

^  ^  i'o  V  4;r(/+m)!  ^  ^  ^ 

which  can  be  converted  to  the  geoid  anornaly  SN  (in  m)  by  dividing  by 
g  =  -GM/R-: 

5N  =  ^=  /?|s(fr*Z  A/^^~^^lT(cos0)[A^coswq)+ B^sinwtp]! 

^  l/=2  m=()  V  47t(/+w)!  I 

The  geoid  anomalies  calculated  here  are  referenced  to  a  theoretical  hydrostatic 
sphere  to  remove  the  effect  of  the  earth's  rotation  (Hager,  1984).  Gravity  is  the 
derivative  of  the  gravitational  potential  with  respect  to  radial  distance,  so  the 
radial  gravity  anomaly  is; 


Gravity  gradient  is  the  derivative  of  gravity  with  respect  to  radial  distance,  so 
the  radial  gravity  gradient  anomaly  is: 


72 


Sln  = 


djSgr) 

dr 


^  CM  ly  (MX/+2)jflf*^y  fghm-m)'- 

R  Is  R-^  4^{'+''0! 


iy(cos0)[A^'cosm(p  +  B/"  sin/?2(p] 


Evaluating  at  r  =  /?  and  using  the  spherical  harmonic  coefficients  2-20  from  the 
GEM-L2  model  (Lerch  et  ai,  1982),  the  equations  simplify  to: 

Z  P'/'(cose)[Afcosm(|)-t-  B^sinmcp]! 

\h2m=()  V  47r(/+w)!  I 


=  ^/^^^2^Ff(cos0)[Afcosw(|)+  B/”sinm(p]j 

/?2  1/^2  m=()  V  4K(/+m)!  I 

srn  =  -^lz  (/+1)(/+2)Z  P;”(cos0)[A;"cos/n(p4  B^sinmcp]) 

l/=2  ,«=()  V  47r(/+m)l  I 

It  is  important  to  note  that  the  GEM-L2  coefficients  must  be  multiplied  by 
before  they  are  plugged  into  these  equations  to  be  consistent  with  the  spherical 
harmonic  normalization  used  in  this  study.  The  three  geophysical  control  data 
sets  are  constructed  by  calculating  the  values  of  the  geoid,  gravity  and  gravity 
gradient  anomalies  at  the  36  feature  locations  (Table  3.3). 


Variation-Di.stance  Relationships 

The  different  characteristics  of  the  contructed  geophysical  data  sets  are 
apparent  in  plots  of  absolute  difference  versus  angle  (Figs.  3.16-3.18).  The 
geoid  plot  shows  a  clean  and  fairly  symtnetric  degree  2  pattern,  with  an  angular 
sampling  distance  of  ~  102".  The  gravity  plot  is  a  little  more  dispersed,  with 
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weaker  symmetry  and  an  angular  sampling  distance  of  -  95°  (expansion  to  degree 
2].  Finally,  the  gravity  gradient  plot  shows  even  more  dispersion  and  an  angular 
sampling  distance  of  -  67°  [expansion  to  degree  3).  A  comparison  of  these  plots 
to  those  for  the  mantle  components  clearly  illustrates  the  complexity  of  the 
geochemistry  data.  Even  the  gravity  gradient  data  (dominated  by  short 
wavelength  energy]  appears  to  have  less  small-scale  variation  (larger  angular 
sampling  distance]  than  all  of  the  mantle  component  data  sets. 

Variation  Reduction  by  Filtering 

The  same  circular  filter  technique  outlined  above  is  applied  to  the 
geophysics  data  to  see  its  effect  (Figs.  3.19-3.21).  The  filtered  geoid  data  set 
retains  its  strong  degree  2  signature  (angular  sampling  distance  -  93°],  but  there 
is  a  slight  increase  in  the  dispersion  of  the  data  points.  Like  the  geoid,  the 
filtered  gravity  data  maintains  its  angular  sampling  distance  (~  93°|  and  it  shows 
a  slight  decrease  in  data  dispersion.  The  gravity  gradient  data  is  most  affected  by 
the  filtering  process.  The  data  dispersion  due  to  large  variation  at  small  and 
large  angles  is  reduced.  In  addition,  the  angular  sampling  distance  is  increased  to 
~  77°.  corresponding  to  spherical  harmonic  expansion  to  degree  2. 


EXPANSION  OF  GEOPHYSICAL  AND  GEOCHEMICAL  DATA  SETS 
By  choosing  the  sampling  distances  based  upon  the  inherent  variation- 
distance  relationships  of  the  different  data  sets,  the  problem  of  aliasing  is 
eliminated.  Of  course,  the  location  patterns  that  result  from  spherical  harmonic 
expansions  may  not  represent  the  true  patterns  as  they  exist  in  the  mantle,  but 
without  a  more  extensive  global  data  set,  there  is  no  way  to  better  approximate 
the  true  pattern.  Coefficients  will  be  found  for  all  six  geophysical  data  sets 


(filtered  and  unfiliered],  for  the  BMII  and  DMM  data  sets  and  for  the  filtered 
EMI  and  HIMU  data  sets. 

Solving  for  the  spherical  harmonic  coefficients  needed  to  expand  a  given 
function  is  a  linear  inverse  problem.  More  specifically,  the  expansion  of  the 
mantle  components  or  geophysical  measures  is  a  di.screte  linear  inverse  problem, 
since  the  data  are  discrete  observations.  The  terminology  and  symbology  used 
here  to  discuss  inverse  problems  is  that  of  Menke  (1989).  Values  of  the  mantle 
components  or  geophysical  measures  at  the  feature  locations  form  a  vector  of 
data  values  d  (/Vxl],  The  unknown  spherical  harmonic  coefficients  form  a 
vector  of  model  parameters  ni  (Mxl  (.  Relating  the  two  is  the  data  kernel  matrix 
G  [/VxA/J,  composed  of  Legendre  polynomials  (functions  of  colatitude(  combined 
with  sine  and  cosines  (functions  of  longitude).  In  matrix  form  the  equation  is: 
(im  =  d,  or  written  out  more  explicitly: 


V  Fq(cos0())  . 

. .  V  l’^;;(cosG())  cosLtpo 

V  P^(cos0())  sinLtpo 

’a?)' 

* 

d() 

P[j(C0SeA.)  .  . 

.  V  I^^(cos9a/)  cosL(p/v 

V  l^(cos0A/)  sinL(pA'_ 

At 

A. 

.d.A'- 

where  L  is  the  maximum  degree  of  the  expansion,  N  is  the  number  of  data 
observations  and  V  is  the  normalization  factor  mentioned  earlier. 

Least  Squares  Method 

Theory.  If  the  equation  (Jm  =  d  provides  enough  information  to  uniquely 
determine  the  model  parameters  or  the  best  fit  to  the  model  parameters,  then 
solving  for  the  spherical  harmonic  coefficients  from  degrees  0  to  5  is  an  even- 
determined  problem  \N  =  .36,  M  =  .36)  and  solving  for  the  coefficients  from 
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degrees  0  to  <5  is  an  overdetermined  problem  \N  =  36,  M  <36],  For  an 
overdetermined  system  of  equations  Gm  =  d,  with  more  equations  than 
unknowns,  there  is  no  exact  solution.  The  least  squares  method  finds  the  model 

parameters  that  minimize  the  error  between  the  observed  data  and  the  predicted 
data,  ie.  it  minimizes  the  L-^  norm  of  the  prediction  error; 

L2  norm:  || e||2  =  'y  ^  | ej^  ,  where  e,-  =  d^”^ 

1=1 

When  solving  for  the  model  parameters  m  {spherical  harmonic  coefficients],  it  is 
best  to  use  QR  decompositon.  The  normal  equations  G^Gm  =  G^d  lead  to  the 
solution;  =  (G^G)'ki^d,  but  if  G^G  is  ill-conditioned,  then  taking  its 
inverse  leads  to  inaccurate  solutions.  QR  decomposition  is  more  accurate  than 
the  normal  equations  for  ill-conditioned  matrices.  It  decomposes  the  data  kernel 
matrix  G  into  two  matrices;  Q  [orthogonal]  and  R  [upper  triangular]:  QRm  = 
d,  with  solutions;  =  R'^Q  ^d. 

Application.  As  a  test  of  the  viability  of  the  least  squares  method,  the 
spherical  harmonic  coefficients  for  the  EMU  percentage  data  and  the  geoid 
anomaly  data  are  solved  for  in  nested  groupings  from  degrees  0-1  up  to  degrees 
0-5.  As  the  data  is  expanded  out  to  greater  degrees,  the  coefficients  should 
decrease  smoothly.  Table  3.4  shows  how  the  degree  2  coefficients  vary  as  the 
tv/o  data  sets  are  expanded  out  to  progressively  higher  degrees.  Only  the  Ao  and 
A2  coefficients  for  the  geoid  and  the  A2  coefficient  for  EMU  decrease  smoothly 
for  the  degrees  0-2  through  degrees  0-4  expansions.  The  other  coefficients 
either  get  larger  or  oscillate.  When  solving  the  even-determined  system  [degrees 
0-5],  all  of  the  coefficients  experience  a  large  increase  or  decrease,  indicating  a 
verv  unstable  solution. 
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Since  the  geoid  coefficients  are  known,  the  correlations  (by  degree) 
between  the  actual  coefficients  and  the  computed  coefficients  for  the  nested 
groupings  can  be  calculated.  The  correlation  coefficient  for  two  sets  of 

coefficients  and  |A2,B2]  is  given  by  the  ratio  of  covariance  to  variance 

at  each  harmonic  degree  (Richards  and  Hager,  1988): 

X  lAlTAlT  +  BVfBlT] 

m={) 

n  -  — ,  -  ■  - 

V  t  [(Air)-  +  (BI?'F)i  ((A2?')-  +  (B27F] 

m=()  ni=() 

Correlations  with  the  actual  geoid  coefficients  can  only  be  made  at  degrees  2  and 
higher  since  the  actual  degree  0  and  1  coefficients  are  zero.  Correlations  of  the 
actual  geoid  coefficients  to  those  calculated  using  least  squares  are: 


Expansion 

Degree  2 

Correlation  Coefficient  [r^] 

Degree  3  Degree  4 

Degree 

Degrees  0-2 

0.046 

--- 

... 

Degrees  0-3 

0.960 

0.794 

... 

Degrees  0-4 

0.884 

0.469 

0.597 

... 

Degrees  0-5 

-0.219 

-0.105 

-0.035 

0.03 1 

The  expansion  for  degrees  0-3  shows  the  best  correlation,  but  there  is  no 
consistency  from  expansion  to  expansion.  Since  the  least  squares  solutions  do  not 
exhibit  consistent,  stable  behavior,  it  appears  that  the  system  (»ni  =  d  does  not 
provide  enough  information  to  uniquely  determine  the  model  parameters  (or  a 
best  estimate  for  them].  This  indicates  that  the  svstem  is  not  even-  or 
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overdetermined,  but  mixed-determined  (neither  completely  overdetermined  nor 
completely  underdetermined]  and  requires  a  more  sophisticated  method  to  solve 
for  the  coefficients. 


Singular  Value  Decompositon  Method 

Theory.  Singular  value  decompositon,  or  SVD,  is  one  way  to  solve  a 
mixed-determined  problem.  Its  purpose  is  to  partition  the  system  of  equations 
into  an  overdetermined  part  (that  can  be  solved  in  the  least  squares  sense]  and  an 
underdetermined  part  (that  can  be  solved  assuming  some  a  priori  information]. 
For  the  general  equation  Gm  =  d,  it  is  like  a  transformation  to  the  system  G'm' 
=  d',  where  m'  is  composed  of  an  overdetermined  part,  m®  and  an 
underdetermined  part  (Menke,  1989): 


Gni  =  d  — >  G'm'  =  d' 


G"'  d 

[mo'' 

d«' 

0  G"'. 

Lm“'. 

d“. 

SVD  decomposes  the  data  kernel  matrix  G  into  three  matrices:  (»  = 
UAV^.  The  matrix  G  is  an  NxN  matrix  of  orthonormal  (orthogonal  and  of  unit 
length]  eigenvectors  that  span  the  data  space  S(d).  Similarly,  the  matrix  V  is  an 
MxM  matrix  of  orthonorma!  vectors  that  span  the  model  parameter  space  S(m). 
The  matrix  A  is  an  NxM  diagonal  eigenvector  matrix  with  nonnegative  diagonal 
elements  called  singular  values,  arranged  in  order  of  decreasing  size.  Some  of 
the  singular  values  may  be  zero,  making  it  easy  to  partition  the  matrix  into  a 
submatrix  A^.  with  p  nonzero  singular  values,  and  several  zero  matrices: 


A  = 


A/,  0 
0  0 


This  simplifies  the  data  kernel  decomposition  to:  G  =  G/,A^V^ 


where  and  V^  are  the  first  p  columns  of  G  and  respectively. 
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For  the  equations  (ilm  =  cl,  the  solution  is:  =  V^A^'U^^d,  called  the 

natural  solution  (Menke,  1989).  If  the  equation  GM  =  d  is  to  some  degree 
underdetermined,  specifies  the  combinations  of  model  parameters  for  which 

the  equation  does  provide  information;  these  combinations  lie  in  a  subspace  of 
the  model  parameter  space  S^(m).  On  the  other  hand,  if  GM  =  d  is  to  some 
degree  overdetermined,  then  A^  specifies  the  combinations  of  model  parameters 

that  the  product  Gm  is  capable  of  resolving;  these  products  span  a  subspace  of 
the  data  space  S^(d).  If  none  of  the  singular  values  are  zero,  there  are 

undoubtedly  some  very  close  to  zero  that  are  affecting  the  solution  variance. 

One  way  to  reduce  the  solution  variance  is  to  select  a  cutoff  size  for  the  singular 
values  and  exclude  any  singular  values  smaller  than  this  lie.  artificially  decide  the 
size  of  p,  the  number  of  nonzero  singular  values],  i  nis  is  equivalent  to  throwing 
away  some  combinations  of  the  model  parameters  (thus  reducing  the  sizes  of 
and  V^j.  However,  if  the  singular  values  excluded  are  small,  then  the  solution 

will  be  clo.se  to  the  natural  solution,  though  the  data  and  model  resolution  will  be 
worse.  This  is  a  classic  trade  off  situation  between  resolution  and  variance 
(Menke,  1989). 

It  is  also  possible  to  dampen  the  smaller  singular  values  instead  of 
throwing  them  away  [equivalent  to  the  damped  least  squares  method].  The 
drawbacks  to  this  method  are  that  the  solution  is  no  longer  close  to  the  natural 
solution,  the  data  and  model  resolution  are  worse  and  the  damping  parameter 
must  be  determined  by  trial-and-error.  For  this  study,  various  methods  are  used 
to  try  to  determine  the  optimum  number  of  singular  values  to  keep  \p]  and  all 
singular  values  with  index  >  p  are  dropped. 

Desired  number  of  singular  values.  The  first  step  in  determining  the 
desired  number  of  singular  values  is  to  look  at  the  data  kernel  spectrums  [plots  of 
the  size  of  the  singular  values  versus  their  index]  for  the  mantle  component  data 
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kernel  and  the  geophysical  data  kernels  (Figs.  3.22-3.25).  For  the  mantle 
components,  »he  data  kernel  (i  is  only  a  function  of  location,  so  it  is  the  same  for 
all  four  components.  For  the  geophysical  data,  the  data  kernels  are  constructed 
differently,  so  that  all  three  equations  with  geoid,  gravity  and  gravity  gradient 
data  are  solving  for  the  same  spherical  harmonic  coefficients.  With  respect  to 
the  mantle  component  data  kernel,  terms  in  the  geoid,  gravity  and  gravity 


GM(/+i) 

gradient  data  kernels  are  multiplied  by  the  additional  factors  of  R,  /^2  and 

-GM(/+i)(/+2) 

,  lespectively. 

For  comparison,  spectrums  for  the  degrees  0-1, 0-2,  0-3,  0-4  and  0-5 
expansions  are  all  plotted,  but  the  emphasis  here  will  be  on  getting  reasonable 
results  using  the  degrees  0-5  expansion.  All  three  geophysical  spectrums  and  the 
geochemical  spectrum  for  this  expansion  show  the  singular  values  gradually 
decreasing  in  value,  with  the  last  five  or  so  singular  values  being  very  close  to 
zero.  There  is  no  obvious  cutoff  size  for  the  singular  values  apparent  in  these 
plots,  so  other  methods  .'iust  be  used  to  estimate  p. 

For  the  geophysical  control  set,  it  is  possible  to  find  the  number  of 
singular  values  p  needed  to  most  closely  approximate  the  actual  coefficients.  The 
root  mean  square  error  between  the  actual  and  estimated  geophysical  coefficients 
is  given  by: 


coefficient  mis  error  = 


X  (mf  ‘  -  nif  >)' 


M 


where  M  is  the  number  of  coefficients  (model  parameters).  A  plot  of  coefficient 
rms  error  versus  the  number  of  singular  values  retained  (Fig.  3.26)  indicates  that 
30.  26  and  14  singular  values  should  be  retained,  for  geoid,  grtivity  and  gntvity 
gradient,  respectively,  to  most  closely  approximate  the  actual  coefficients.  'I'hcse 
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values  are  indicated  on  the  data  kernel  spectrum  plots  (Figs.  3.22-3.24).  It  is 
important  to  note  that  the  more  a  field  is  dominated  by  high  degree  energy,  the 
fewer  singular  values  it  takes  for  the  rms  error  to  explode  [at  least  for  these 
sparse  data  sets]. 

Since  the  coefficients  for  the  geochemistry  data  are  not  known,  there  is  no 
way  to  measure  how  closely  the  estimated  coefficients  match  the  actual 
coefficients.  What  can  be  done  is  to  try  to  match  the  observed  data  as  closely  as 
possible,  while  keeping  the  solution  variance  at  a  minimum.  As  a  first  step, 
trade-off  curves  are  constructed  to  bracket  the  range  of  p  values  that  balance  the 
size  of  the  model  variance  and  the  spread  of  the  model  resolution  (Figs.  3.27- 
3.30).  The  size  of  the  model  variance  is  based  upon  the  unit  covariance  matrix 
of  the  model  parameters,  which  characterizes  the  degree  of  error  amplification 
that  occurs  in  the  mappinp  from  data  to  model  parameters  (Menke,  1989). 
Assuming  that  the  data  within  the  four  mantle  component  vectors  and  the  three 
geophysical  vectors  are  uncorrelated  and  have  uniform  variance  C5^  [a  reasonable 
assumption  for  the  mantle  component  vectors  based  upon  the  findings  in  Chapter 
2],  the  covariance  matrix  of  the  model  parameters  is  given  by; 

[cov  mcst]  =  G'®[cov  d]G'^^= 

where  G'^  is  the  generalized  inverse,  which  for  singular  value  decomposition  is: 


G-s  =  V^A;,iuJ 


The  unit  covariance  matrix  is: 

[covu  mcst]  =  mc-si]  =  G-eG-g"^  = 

Finally,  the  size  of  the  model  variance  is; 

.  _ ,,2  M  M 

sizB([cov  m“‘])  =  II V varu  lb  =  ^  [vary  ^  [covy 

1=1  1=1 


where  M  is  the  number  of  model  parameters.  To  summarize,  the  size  of  the 
model  variance  is  the  sum  of  the  variances  of  the  model  parameters,  which  arc 
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the  diagonal  elements  of  the  model  parameter  covariance  matrix.  With 
increasing  values  of  p,  the  size  of  the  model  variance  will  increase. 

Since  resolution  is  optimal  when  the  resolution  matrices  are  identity 
matrices,  it  is  possible  to  quantify  the  spread  of  model  resolution  based  on  the 
size  of  the  off-diagonal  elements  of  the  model  resolution  matrix  R  (Menke, 
1989): 

M  M 

spread(R)-||R-llg  =  XI[Ro-J// 

/=!  j=\ 

where  I  is  the  identity  matrix  and  R  =  ,  [m^si-  Rm^*^].  With  increasing 

values  of  p,  the  spread  of  the  model  resolution  will  decrease. 

Trade-off  curves  of  size  of  model  variance  versus  spread  of  model 
resolution,  as  a  function  of  the  number  of  singular  values  retained,  show  two 
asymptotes  (retaining  all  36  singular  values  gives  the  largest  model  variance  size] 
(Figs.  3.27-3.30).  The  ideal  range  for  p,  to  balance  the  two  measures,  is  in  the 
transition  between  the  asymptotes  (Table  3.5). 

Another  way  to  try  and  pin  down  the  desired  number  of  singular  values 
(to  most  closely  approximate  the  data]  is  to  look  at  plots  of  model  rms  error  and 
a  variance  measure  versus  the  number  of  singular  values  retained  (Figs.  3.31- 
3.34).  Model  rms  error  is  given  by: 

model  rms  error  = 

where 

While  VjV^  and  arc  the  identity  matrix,  and  U^LlJ  are  not 

necessarily  the  identity  matrix,  since  and  do  not  in  general  span  the 

complete  data  and  model  spaces  (Menke,  1989).  I'he  variance  measure  used  is: 
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variance  me^isure  =  ^  [Ap^]„- 
/=1 

since  the  solution  variance  is  proportional  to  Ap  .  Again  ,  the  goal  is  to  use  the 
plots  of  these  two  quantities  to  select  p  so  that  the  model  rms  error  and  the 
solution  variance  are  balanced  (Table  3.5). 

Choosing  ranges  for  p  using  trade-off  curves  and  the  model  rms/variance 
curves  is  a  subjective  process.  The  ranges  of  values  are  chosen  by  eye  and  there 
is  no  objective  way  to  select  an  optimal  value  of  p  from  these  methods.  To  make 
the  process  more  objective,  Jacobson  and  Shaw  (1991)  suggest  applying  a 
sequential  F-test  to  SVD  problems  to  find  the  statistically  optimal  solution. 

Given  a  null  model  with  q  parameters  and  a  larger  general  model  with  h 
parameters  \h  >  q\,  testing  the  null  hypothesis  that  the  additional  [h  -  q] 
parameters  in  the  general  model  do  not  improve  the  fit  to  the  data  (compared  to 
the  null  model]  requires  tfic  use  of  the  F-statistic: 

_  ( RSS^  -  RSS/^)  (/2  -  h) 

where  RSS^  and  RSS/^  denote  the  residual  sum  of  squares  for  the  null  and  general 

models,  respectfully,  and  n  is  the  total  number  of  parameters.  F  has  an  F- 

distribution  with  {b  -  q,n  -  b)  degrees  of  freedom.  The  residua!  sum  of  squares 

for  a  given  model  is  defined  as: 

RSS  =  X  (d?  - 
/=! 

Values  of  F  can  be  converted  into  the  probability  that  the  null  hypothesis  is  true, 
ic.  that  the  extra  parameters  do  not  result  in  a  better  fit.  Then  the  quantity  [1  - 
prob(null  hypothesis  true)]  is  the  significance  level  of  the  additional  parameters. 
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For  SVD,  the  sequential  F-test  starts  by  testing  the  significance  of  a  model 
retaining  one  singular  value  against  a  model  retaining  no  singular  values,  then 
continues  to  test  models  retaining  incrementally  more  singular  values  against  the 
current  null  model.  When  a  model  has  reached  tlie  95%  significance  level 
[chosen  for  this  application]  or  higher,  it  becomes  the  null  model  against  which 
subsequent  models  are  to  be  tested,  until  another  model  also  reaches  or  surpasses 
95%  significance  and  takes  its  place.  Figures  3.35-3.41  show  the  F-test  results 
for  the  geophysical  and  geochemical  data  sets  and  Table  3.5  lists  the  resulting 
optimal  p  values.  In  general,  it  appears  that  the  smoother  functions  (longer 
wavelength]  have  higher  numbers  of  significant  singular  values. 

For  determining  the  value  of  />,  the  three  different  methods  agree  quite 
well  (Table  3.5).  The  trade-off  curves  define  the  largest  interval  for  p,  which  is 
constrained  further  by  the  model  rms/variance  curves.  For  every  data  set, 
except  Altered  gravity,  the  value  of  p  determined  by  the  F-test  falls  within  the 
cho.sen  range  of  the  model  rms/variance  curves.  Even  so,  the  F-test  p  value  for 
fiiterexl  gravity  does  not  fall  far  outside  the  model  rms/variance  range  Ip  =  29 
compared  to  25]  and  it  does  fall  within  the  trade-off  range.  Since  the  F-test  p 
values  are  in  agreement  with  the  other  methods  and  are  by  far  the  most  objective 
estimate  from  the  three  methods,  these  values  will  be  used  in  calculating  the 
spherical  harmonic  coefficients. 

Appiication.  How  well  the  estimated  spherical  harmonic  coefficients  of 
the  constructed  geophysical  data  .sets  correlate  with  the  actual  GEM-L2 
coefficients  is  an  indicator  of  how  closely  the  estimated  geochemistry  coefficients 
may  be  expected  to  approximate  their  true  coefficients.  Three  sets  of 
geophysical  SVD  coefficient  solutions  are  all  correlated  with  the  GEM-L2 
coefficients:  those  that  minimize  the  coefficient  rms  error  and  those  that 
minimize  the  model  error  [selected  p  values  from  the  F-test]  for  the  filtered  and 
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unfiitcred  data  sets  ('Fable  3.6).  Remember  that  the  data  kernel  matrices  (»  for 
gravity  and  gravity  gradient  are  modified  so  that  their  spherical  harmonic 
coefficients  are  also  estimates  of  the  GEM-L2  coefficients.  The  correlation 
coefficients  are  calculated  as  outlined  above.  Plots  of  versus  degree  include 

confidence  levels  based  upon  a  student’s  r-test.  The  test  statistic  for  the  /-test  is: 

rj,  _  iiin-  2  _  Tiill 

where  n  is  the  number  of  coefficients  at  that  particular  degree  \{n  -  2)  =  2/J.  T 
has  a  /-distribution  w'ith  {u  -  2)  or  21  degrees  of  freedom.  Given  a  desired 

significance  level  and  the  degrees  of  freedom,  the  value  of  T  can  be  looked  up  in 
a  table.  Then  the  vaiue  that  Vj  should  have  to  achieve  that  significance  level  can 

be  calculated  and  plotted  as  confidence  levels: 


For  the  plots  of  r^  versus  /,  the  geophysical  coefficients  estimated  by 
minimizing  the  coefficient  rms  error  correlate  better  than  those  estimated  by 
minimizing  the  model  error  and,  of  those,  the  unfiltered  data  set  co’.relatcs  better 
than  the  filtered  data  set.  All  three  .sets  of  coefficients  correlate  well  with  the 
actual  GHM-T2  coefficients  at  degree  2,  except  for  filtered  gravity  (F'igs.  3,42- 
3.44  ).  In  all  cases,  the  geoui  coefficient  estimates  correlate  the  best.  In  general, 
gravity  atid  gravity  gradient  correlate  better  at  even  degrees,  with  the  exception 
(4  the  filtered  cociTicients.  For  the  mantle  c<’)mponeni  coefficients,  all  this 
implies  that  the  degree  2  coefficients  arc  probably  good,  but  beyond  that  there  is 
no  guarantee.  Of  tltc  lour  mantle  component  percentage  data  sets  that  are 
expanded,  the  .iltcred  llI.Mli  data  set  is  unique  in  that  it  most  closely  resembles 
the  geoid  data  set  in  the  vari.nion-distancc  plots  (F'igs.  3.13  and  3.16).  Thus, 
there  is  a  goo<1  possii  i!it\  that  at  Ie;!si  the  degree  3  coefficients  for  this  data  set 
are  rea^on:  ble  as  well. 
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Correlation  coefficients  for  the  actual  GEM-L2  coefficients  and  the 
estimated  coefficients  cannot  be  calculated  at  degrees  0  and  1  because  those 
GEM-L2  coefficients  are  equal  to  zero.  In  contrast,  the  estimates  of  these 
coefficients  from  the  constructed  geophysical  data  sets  are  all  positive  numbers 
the  same  order  of  magnitude  as  the  rest  of  the  estimated  coefficients.  This 
di.screpancy  is  caused  by  a  sampling  bias  due  to  the  fact  that  the  oceanic  islands 
are  all  hotspot  related  and  hotspots  are  associated  with  geoid  highs  (Richards  et 
al.,  1988];  no  geoid  lows  are  sampled  to  balance  these  highs.  It  is  unclear  how 
this  bias  may  affect  the  estimates  of  the  other  coefficients. 

The  continuous  layer  model  degree  2  "functions"  for  the  constructed  geoid 
data  set  and  the  mantle  component  percentages  are  reconstructed  on  a  five  degree 
grid  over  the  globe  from  1O<0<17O  and  -18()<(p<18()  using  the  calculated 
coefficients  and  the  appropriate  equations  (Figs.  3.45-3.49).  It  should  be  noted 
that  the  contoured  values  arc  not  actual  geoid  anomaly  values  or  component 
percentages,  but  are  deviations  from  the  average  (degree  0|  geoid  anomaly  value 
or  component  percentage  (average  constructed  geoid  =  13.7  m;  average  filtered 
EMI  =  0.27;  average  EMU  =  0.17;  average  filtered  HIMU  =  0.31;  average  DMM 
=  0.25].  For  comparison,  the  actual  degree  2  geoid  is  constructed  in  the  same 
way  using  the  GEM-L2  coefficients  (average  geoid  =  0.0  m]  (Fig.  3.50).  The 
constructed  geoid  field  agrees  well  with  the  actual  degree  2  geoid,  as  already 
indicated  by  the  correlation  coefficients.  For  the  mantle  components,  HIMU 
resembles  the  actual  geoid  field  with  two  essentiallv  equatorial  highs  in 
approximately  the  same  locations;  EMI  and  EMU  also  have  two  highs  that 
undulate  above  and  below  the  equator  with  a  longitudinal  shift  of  -35'’  to  the  east 
with  respect  to  the  actual  geoid  ]EMI1  has  less  offset  than  EMI];  and  DMM.  with 
its  two  highs  and  two  lows  resembles  none  of  the  other  degree  2  expansions 
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Of  all  the  mantle  component  data  sets,  filtered  HIMU  has  the  best  chance 
of  gettinji  reasonable  valecs  for  the  degree  3  coefficients.  I'hc  2-3 

function  for  filtered  HIMU  is  reconstructed  as  l:>efore  (Fig.  3.51).  This  can  be 
compared  to  the  degrees  2-3  geoid  reconstructed  from  the  GEM-L2  coefficients 
(Fig.  3.52). 


SUMMARY 

Viewing  the  distribution  of  the  OIB  reservoir  as  a  continuous  layer  in  the 
mantle  and  using  approximation  methods  to  solve  for  the  spherical  harmonic 
coefficients  of  its  expansion  reveals  the  following: 

•  The  mantle  end-member  component  percentage  data  have  a  lot  of  short 
wavelength  energy  relative  to  equally  limited  geoid,  gravity  and 
gravity  gradient  control  data  sets. 

•  With  the  currently  available  data,  solving  for  the  spherical  htirmonic 
coefficients  is  a  mixed-determined  problem,  requiring  the  use  of 
singular  value  decomposition  (SVD|  to  get  viable  solutions. 

•  'I'he  F-test  is  a  simple,  objective  way  to  determine  the  number  of 
singular  values  to  retain  in  SVD  for  the  statistically  optimal  solution. 

•  With  the  current  data  coverage,  only  the  degree  2  spherical  harmonic 
coefficients  can  be  estimated  with  a  reasonable  level  of  confidence 
using  SVD. 


Continuous  layer  model  degree  2  HIMU  closely  resembles  the  degree  2 
geoid. 

Continuous  layer  model  degree  2  CMI  and  EMU  resemble  a  longitude- 
shifted,  undulating  degree  2  geoid. 

Continuous  layer  model  degree  2  DMM  does  not  resemble  the  degree  2 
geoid  or  the  degree  2  expansion  of  any  other  mantle  component. 
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Table  3.1.  Mantle  end-member  component  percentage.s*  for  the  average  isotopic  signatures  of  the  geographic  features 
I  island  groups,  islands,  ridges,  seamounts]  represented  m  the  OIB  data  set  with  their  locations  and  the  number  of  samples 
for  each  feature  [in  braces]. 
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Table  3.1.  Continued. 
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Table  3.2.  Separation  of  the  OIB  feature  data  set  into  27  Dupal-type  features  and 
9  DMM-type  features,  based  upon  the  percentage  of  the  DMM  mantle 
component. 


Dupal-type  Features 

%DMM 

Gough 

3.20 

Shimada  Seamount 

6.01 

Tristan  de  Cunha 

6.70 

San  Felix/San  Ambrosio 

8.08 

Walvis  Ridge 

1  l.SO 

Azores 

13.74 

d'rinidade 

14.84 

Cook-Austral  Islands 

16.22 

Tubuai-Austral  Islands 

16.93 

New  England  Seamounts 

17.00 

Christmas 

\1.25 

St.  Helena 

17.2.S 

Kerguelen  Plateau 

17.39 

Samoa  Ishinds 

18.88 

Comores  Archipelago 

19.37 

Fernando  de  Noronha 

19.80 

Marquesas  Archipelago 

20.72 

Pitcairn 

21..S7 

Cameroon  Line 

21.7.^ 

Sala  Y  Gomez 

23.18 

Society  Ridge 

23.31 

Cape  Verde  Islands 

26.06 

Juan  Fernandez  Islands 

26.63 

Crozet  Islands 

27.34 

Mascareignes 

28.87 

Amsterdam/St.  Paul 

29.20 

Louisville  Seamount  Chain 

3 1 .84 
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1'able  3.2.  Continued. 


DMM-type  Features 

%DMM 

Balleny 

32.17 

Cocos 

38.53 

Marion/Prince  Edward 

40.88 

Galapagos  Islands 

41.32 

Ascension 

41.69 

Ponape 

46.14 

Hawaiian  Islands 

47.18 

Iceland 

52.81 

Nunivak 

58.73 

Table  3.3.  GcimcI,  gravity  and  gravity  gradient  tinonialy  values,  calculated  using  the  degrees  2-20  Gli.M-L2  spherictil 
harmonic  coctYicienis,  at  the  geographic  features  [island  groups,  ishinds,  ridges,  seamounts]  represented  in  the  OIB  dai;t 
set  with  their  locations. 
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Table  3.3.  Continued. 
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Table  3.4.  Change  in  degree  2  spherical  harmonic  coefficients  for  the  EMU 
percentage  data  and  the  geoid  anomaly  data  as  the  data  sets  are  expanded  to 
progressively  higher  degrees. 


Expansions 

A2 

A2 

Bi 

A  2 

A2 

Geoid 

Degrees  0-2 

-1.491E-05 

-1.049E-06 

-4.116E06 

l.OOlE-05 

-6.258E-06 

Degrees  0-3 

-1.511E-05 

2.881E-06 

1.551E-06 

6.610E-06 

-8.569E-06 

Degrees  0-4 

-1.825E-05 

-7.480E-06 

3.639E-06 

6.044E-O6 

-1.108E-06 

Degrees  0-5 

3.946E-05 

3.357E-05 

-4.237E-04 

-1.436E-04 

-3.682E-06 

EMU 

Degrees  0-2 

-0.065243 

-0.063197 

0.269622 

0  108911 

-0.008334 

Degrees  0-3 

-0.076403 

-0.203182 

0.466347 

0.148947 

0.053770 

Degrees  0-4 

-0.170987 

-0.221924 

0.785142 

0.167714 

0.087842 

Degrees  0-5 

2.543222 

10.864696 

-30.064709 

-20.805890 

-10.644588 

Table  3.5.  Optimal  values  or  ranges  of  values  for  p  [the  number  of  singular  values  retained)  for  the  best  approximations 
of  t..e  observed  data  that  keep  solution  variance  to  a  minimum,  as  determined  by  three  ditfereni  methods:  trade-oil 
curvis.  model  rms  error  and  variance  curves,  and  the  F-test. 
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Table  3.6.  Summary  of  correlation  coefficients  between  the  GhM-L2 
coefficients  and  three  sets  of  estimated  geophysical  coefficients  that  minimize  the 
coefficient  rms  error  and  that  minimize  the  model  error  for  filtered  and 
unfiltered  data  sets. 


Data  Set 

p  Value 

Degree  2 

Degree  3 

Degree  4 

Degree  5 

Minimizing  coefficient  ] 

rms  error 

geoid 

30 

0.988 

0.909 

0.744 

0.705 

gravity 

26 

0.981 

0.605 

0.734 

0.310 

gravity  gradient 

14 

0.926 

0.399 

0.457 

-0.014 

Minimizing  model  error  (F-test)  - 

unfiltered 

geoid 

25 

0.988 

0.863 

0.726 

-0.031 

gravity 

20 

0.891 

0.264 

0.394 

-0.145 

gravity  gradient 

20 

0.779 

-0.084 

0.543 

0.085 

Minimizing  model  error  (F-test)  - 

filtered 

geoid 

24 

0.629 

0.537 

-0.256 

0.016 

gravity 

29 

0.356 

0.723 

-0.331 

-0.333 

gravity  gradient 

24 

0.709 

0.334 

-0.720 

-0.491 
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Fig.  3.1.  Variation-distance  plot  for  the  EMI  mantle  component  showing  the 
range  of  variation  in  the  component  percentage  with  angular  distance  between 
the  feature  locations.  To  account  for  the  variation  requires  a  minimum  sampling 
distance  of  ~  14.5°  [degree  12  expansion]. 
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Fig.  3.2.  Variation-distance  plot  for  the  EMU  mantle  component  showing  the 
range  of  variation  in  the  component  percentage  with  angular  distance  between 
the  feature  locations.  To  account  for  the  variation  requires  a  minimum  sampling 
distance  of  ~  39°  [degree  4  expansion]. 
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Fig.  3.3.  Variation-distance  plot  for  the  HIMU  mantle  component  showing  the 
range  of  variation  in  the  component  percentage  with  angular  distance  between 
the  feature  locations.  To  account  for  the  variation  requires  a  minimum  sampling 
distance  of  -  14.5°  [degree  12  expansion]. 
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Fig.  3.4  Variation-distance  plot  for  the  DMM  mantle  component  showing  the 
range  of  variation  in  the  component  percentage  with  angular  distance  between 
the  feature  locations.  To  account  for  the  variation  requires  a  minimum  sampling 
distance  of  ~  57°  [degree  3  expansion]. 


Fig.  3.5.  Three-dimensional  plot  of  the  geographic  feature  principal  component 
data.  Axes:  X  =  Zl,  Y  =  Z2,  Z  =  Z3.  Symbols:  -i-  =  features  inside  the 
DUPAL  belt,  o  =  features  outside  the  DUPAL  belt. 
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Fig.  3.6.  Three-dimensional  plot  of  the  geographic  feature  principal  component 
data,  with  symbols  distinguishing  percentages  of  the  DMM  component.  Axes:  X 
=  Zl,  Y  =  Z2,  Z  =  Z3.  Symbols;  -t-  =  <10%  DMM,  x  =  10-20%  DMM, 
diamond  =  20-30%  DMM,  square  =  30-40%  DMM,  o  =  40-50%  DMM,  A  = 
>50%  DMM.  Most  striking  is  the  large  spatial  separation  in  the  30-40% 
category  between:  (1 )  Louisville  =  31.84%,  12)  Dalleny  =  32.17%,  [3]  Cocos  = 
38.5V^. 
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I  ig.  3.7.  Variation-distance  plot  for  the  EMI  mantle  component  for  the  DUPAL 
features  only  [<32%  DMM],  showing  the  range  of  variation  in  the  component 
{percentage  with  angular  distance  between  the  feature  locations.  Using  the 

DUPAL  features  only  shows  no  reduction  in  the  small-scale  variation  for  the 
EM!  component. 
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Fig.  3.8.  Variation-distance  plot  for  the  EMU  mantle  component  for  the 
DUPAL  features  only  (<32%  DMM],  showing  the  range  of  variation  in  the 
component  percentage  with  angular  distance  between  the  feature  locations. 
Using  the  DUPAL  features  only  shows  no  reduction  in  the  small-scale  variation 
for  the  EMU  component. 
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Fig.  3.9.  Variation-distance  plot  for  the  HIMU  mantle  component  for  the 
DUPAL  features  only  [<32%  DMM],  showing  the  range  of  variation  in  the 
component  percentage  with  angular  distance  between  the  feature  locations. 
Using  the  DUPAL  features  only  shows  no  reduction  in  the  small-scale  variation 
for  the  HIMU  component. 
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Fig.  3.10.  Variation-distance  plot  for  the  DMM  mantle  component  for  the 
DUPAL  features  only  [<32%  DMM],  showing  the  range  of  variation  in  the 
component  percentage  with  angular  distance  between  the  feature  locations. 
Using  the  DUPAL  features  only  does  show  a  reduction  in  the  small-scale 
variation  for  the  DMM  component,  with  an  increase  in  minimum  samplinc 
distance  from  ~  57°  to  89°. 
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Fig.  3.1 1.  Variation-distance  plot  for  the  filtered  EMI  data  set,  showing  the 
range  of  variation  in  the  component  percentage  with  angular  distance  after  the 
circular  filter  is  applied.  The  result  is  a  reduction  in  the  small-scale  variation, 
with  an  increase  in  minimum  sampling  distance  from  ~  14.5°  to  37°  [expansions 
to  degrees  12  and  5,  respectively]. 
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Filtered  EMU 


Fig.  3.12.  Variation-distance  plot  for  the  filtered  EMU  data  set,  showing  ihe 
range  of  variation  in  the  component  percentage  with  angular  distance  after  the 
circular  filter  is  applied.  The  result  is  an  increase  in  the  small-scale  variation, 
with  a  decrease  in  minimum  sampling  distance  from  ~  39°  to  27°  [expansions  to 
degrees  4  and  7,  respectively]. 


Filtered  HIMU 


Fig.  3.13.  Variation-distance  plot  for  the  filtered  HIMU  data  set,  showing  the 
range  of  variation  in  the  component  percentage  with  angular  distance  after  the 
circular  filter  is  applied.  The  result  is  a  dramatic  decrease  in  the  small-scale 
variation,  with  an  increase  in  minimum  sampling  distance  from  ~  14.5°  to  83° 
[expansions  to  degrees  12  and  2,  respectively]. 


niicred  DMM 


Fig.  3.14.  Variation-distance  plot  for  the  filtered  DMM  data  set,  showing  the 
range  of  variation  in  the  component  percentage  with  angular  distance  after  the 
circular  filter  is  applied.  The  result  is  a  decrease  in  the  small-scale  variation, 
with  an  increase  in  minimum  sampling  distance  from  ~  57°  to  102°  [expansions 
to  degrees  3  and  2,  respectively]. 
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Fig.  3.16.  Variation-distance  plot  for  the  constructed  geoid  data  set  showing  the 
range  of  variation  in  the  geoid  with  angular  distance  between  the  feature 
locations.  To  account  for  the  variation  requires  a  minimum  sampling  distance  of 
~  102°  [degree  2  expansion]. 
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Fig.  3.17.  Variation-distance  plot  for  the  constructed  gravity  data  set  showing 
the  range  of  variation  in  gravity  with  angular  distance  between  the  feature 
locations.  To  account  for  the  variation  requires  a  minimum  sampling  distance  of 
~  95°  [degree  2  expansion]. 
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Fig.  3.18.  Variation'distance  plot  for  the  constructed  gravity  gradient  data  set 
showing  the  range  of  variation  in  the  gravity  gradient  with  angular  distance 
between  the  feature  locations.  To  account  for  the  variation  requires  a  minimum 
sampling  distance  of  -  67°  (degree  3  expansion). 
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Fig.  3.19.  Variation-distance  plot  for  the  filtered  geoid  data  set,  showing  the 
range  of  variation  in  the  geoid  with  angular  distance  after  the  circular  filter  is 
applied.  The  filtered  geoid  data  set  retains  essentially  the  sa’ne  angular  sampling 
distance  93°]. 
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f  iltered  Gravity 


Fig.  3.20.  Variation-distance  plot  for  the  filtered  gravity  data  set,  showing  the 
range  of  variation  in  gravity  with  angular  distance  after  the  circular  filter  is 
applied.  The  filtered  gravity  data  set  retains  essentially  the  same  angular 
sampling  distance  [~  93°]. 
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Fig.  3.21.  Variation-distance  plot  for  the  filtered  gravity  gradient  data  set, 
showing  the  range  of  variation  in  the  gravity  gradient  with  angular  distance  after 
the  circular  filter  is  applied.  Filtering  reduces  the  small-scale  varition  in  the 
gravity  gradient  data  set,  with  an  increase  of  angular  sampling  distance  from  ~ 
67°  to  77°  [expansions  to  degrees  3  and  2,  respectively). 
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Geoid  Data  Kernel  Spectrum 
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Fig.  3.22.  Data  kernel  spectrums  for  the  constructed  geoid  data  kernel  G. 
Symbols  for  the  different  expansions:  •  =  degrees  0-1,  +  =  degrees  0-2,  *  = 
degrees  0  3,  o  =  degrees  0-4,  x  =  degrees  0-5.  For  the  degrees  0-5  expansion, 
the  singular  values  approach  zero,  but  there  is  no  obvious  cutoff  value. 
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Fig.  3.23.  Data  kernel  spectrums  for  the  constructed  gravity  data  kernel  G. 
Symbols  for  the  different  expansions:  •  =  degrees  0-1,  +  =  degrees  0-2,  *  = 
degrees  0-3,  o  =  degrees  0-4,  x  =  degrees  0-5.  For  the  degrees  0-5  expansion, 
the  singular  values  approach  zero,  but  there  is  no  obvious  cutoff  value. 


40 


Fig.  3.24.  Data  kernel  spectrums  for  the  constructed  gravity  gradient  data 
kernel  G.  Symbols  for  the  different  expansions:  •  =  degrees  0-1,  +  =  degrees  0- 
2,  *  =  degrees  0-3,  o  =  degrees  0-4,  x  =  degrees  0-5,  For  the  degrees  0-5 
expansion,  the  singular  values  approach  zero,  but  there  is  no  obvious  cutoff 
value. 
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Fig.  3.25.  Data  kernel  spectrums  for  the  mantle  component  data  kernel  G. 
Symbols  for  the  different  expansions;  •  =  degrees  0-1,  +  =  degrees  0-2,  *  = 
degrees  0-3,  o  =  degrees  0-4,  x  =  degrees  0-5.  For  the  degrees  0-5  expansion, 
the  singular  values  approach  zero,  but  there  is  no  obvious  cutoff  value. 
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Fig.  3.26.  Plot  of  the  root  mean  square  error  [rms  error],  as  a  function  of  the 
number  of  singular  values  retained,  between  the  actual  GEM-L2  coefficients  and 
those  coefficients  estimated  by  the  constructed  geoid,  gravity  and  gravity 

gradient  data  sets.  Line  symbols:  -  .  -  .  =  geoid, - =  gravity,  •  ■  •  •  = 

gravity  gradient, - =  root  mean  square  of  the  GEM-L2  coefficients.  P  values 

minimizing  coefficient  rms  error:  geoid  =  30,  gravity  =  26,  gravity  gradient  = 
14. 
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Fig.  3.27.  Trade-off  curve  between  model  variance  and  model  resolution,  as  a 
function  of  the  number  of  singular  values  retained,  for  the  constructed  geoid  data 
set.  Range  for  p  that  balances  the  two  measures  is:  15^<30.  Note  that  trade¬ 
off  curves  are  determined  by  the  data  kernel  matrices  and  so  are  the  same  for 
filtered  and  unfiltered  data  sets. 
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Fig.  3.28.  Trade-off  curve  between  model  variance  and  model  resolution,  as  a 
function  of  the  number  of  singular  values  retained,  for  the  constructed  gravity 
dataset.  Range  for /?  that  balances  the  two  measures  is:  9^<29.  Note  that 
trade-off  curves  are  determined  by  the  data  kernel  matrices  and  so  are  the  same 
for  filtered  and  unfiltered  data  sets. 
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Fig.  3.29.  Trade-off  curve  between  model  variance  and  model  resolution,  as  a 
function  of  the  number  of  singular  values  retained,  for  the  constructed  gravity 
gradient  data  set.  Range  for  p  that  balances  the  two  measures  is:  8</7<26.  Note 
that  trade-off  curves  are  determined  by  the  data  kernel  matrices  and  so  are  the 
same  for  filtered  and  unfiltered  data  sets. 
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Fig.  3.30.  Trade-off  curve  between  model  variance  and  model  resolution,  as  a 
function  of  the  number  of  singular  values  retained,  for  the  mantle  component 
data  set.  Range  for p  that  balances  the  two  measures  is;  I5<p^0.  Note  that 
trade-off  curves  are  determined  by  the  data  kernel  matrices  and  so  are  the  same 
for  filtered  and  unfiltered  data  sets. 
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Fig.  3.31.  Plot  of  model  root  mean  square  error  [rms  error],  as  a  function  of 
the  number  of  singular  values  retained,  between  the  observed  geoid  data  and  the 
geoid  data  predicted  from  the  calculated  coefficients.  Balancing  the  model  rms 
error  and  the  model  variance  gives  this  range  of  p:  2\<p<25  (filtered  and 

unfiltered).  Line  symbols:  - =  unfiltered  model  rms  error, - =  filtered 

model  rms  error,  o — o —  =  model  variance. 
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Fig.  3.32.  Plot  of  model  root  mean  square  error  [rms  error],  as  a  function  of 
the  number  of  singular  values  retained,  between  the  observed  gravity  data  and 
the  gravity  data  predicted  from  the  calculated  coefficients.  Balancing  the  model 
rms  error  and  the  model  variance  gives  this  range  of  p:  20<p<25  (filtered  and 

unfiltered).  Line  symbols:  - =  unfiltered  model  rms  error, - =  filtered 

model  rms  error,  o — o —  =  model  variance. 
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Fig.  3.33.  Plot  of  model  root  mean  square  error  [rms  error],  as  a  function  of 
the  number  of  singular  values  retained,  between  the  observed  gravity  gradient 
data  and  the  gravity  gradient  data  predicted  from  the  calculated  coefficients. 
Balancing  the  model  rms  error  and  the  model  variance  gives  this  range  for  p\ 

20<p<25  (filtered  and  unfiltered).  Line  symbols;  - =  unfiltered  model  rms 

error, - -  filterfid  model  rms  error,  o — o —  =  model  variance. 
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Fig.  3.34.  Plot  of  model  root  mean  square  error  [rms  error],  as  a  function  of 
the  number  of  singular  values  retained,  between  the  observed  mantle  component 
data  and  the  mantle  component  data  predicted  from  the  calculated  coefficients. 
Balancing  the  model  rms  error  and  the  model  variance  gives  this  range  for  p\ 
\6<p<2\  (filtered  EMI),  16^<20  (EMU),  18^<23  (filtered  HIMU),  19<p<22 

(DMM).  Line  symbols;  - =  filtered  EMI, - =  EMU,  •  ■  •  ■  =  filtered 

HIMU,  -  -  =  DMM,  0 — 0 —  =  model  variance. 
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Geoid  and  Filtered  Geoid 


P 


Fig.  3.35.  Plot  of  F-test  significance  level  as  a  function  of  the  number  of 
singular  values  retained  for  the  geoid  and  filtered  geoid  data  sets.  Basically,  the 
test  determines  whether  additional  parameters  [singular  values]  make  a 
significant  contribution  to  the  model  fit  of  the  observed  data  values.  Optimal  p 
values  [for  95%  significance]  are:  p  =  25  [geoid]  and  p  =  24  [filtered  geoid]. 

Line  symbols:  - =  geoid,  -  •  -  ■  =  filtered  geoid, - =  95%  significance 

level. 
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Fig.  3.36.  Plot  of  F-tcsl  significance  level  as  a  function  of  the  number  of 
singular  values  retained  for  the  gravity  and  filtered  gravity  data  sets.  Basically, 
the  test  determines  whether  additional  parameters  [singular  values]  make  a 
significant  contribution  to  the  model  fit  of  the  observed  data  values.  Optimal  p 
values  [for  95%  significance]  are:  p  =  20  [gravity]  and  p  =  29  [filtered  gravity]. 

Line  symbols:  - =  gravity,  -  •  -  =  filtered  gravity, - =  95% 

significance  level. 
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Gravity  Gradient  and  Filtered  Gravity  Gradient 


P 


Fig.  3.37.  Plot  of  F-test  significance  level  as  a  function  of  the  number  of 
singular  values  retained  for  the  gravity  gradient  and  filtered  gravity  gradient 
data  sets.  Basically,  the  test  determines  whether  additional  parameters  [singular 
values]  make  a  significant  contribution  to  the  model  fit  of  the  observed  data 
values.  Optimal  p  values  [for  95%  significance]  are:  p  =  20  [gravity  gradient] 

and  p  =  24  [filtered  gravity  gradient].  Line  symbols:  - =  gravity  gradient, 

-  •  -  ■  =  filtered  gravity  gradient, - =  95%  significance  level. 
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Fig.  3.38.  Plot  of  F-test  significance  level  as  a  function  of  the  number  of 
singular  values  retained  for  the  filtered  EMI  data  set.  Basically,  the  test 
determines  whether  additional  parameters  [singular  values]  make  a  significant 
contribution  to  the  model  fit  of  the  observed  data  values.  For  filtered  EMI,  the 
optimal  p  value  [for  95%  significance]  is:  p  =  20.  Line  symbols:  o — o  =  filtered 
EMI, - =  95%  significance  level. 
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Fig.  3.39.  Plot  of  F-test  significance  level  as  a  function  of  the  number  of 
singular  values  retained  for  the  EMU  data  set.  Basically,  the  test  determines 
whether  additional  parameters  [singular  values]  make  a  significant  contribution 
to  the  model  fit  of  the  observed  data  values.  For  EMU,  the  optimal  p  value  [for 

95%  significance]  is;  p  =  16.  Line  symbols:  o — o  =  EMU, - =  95% 

significance  level. 
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Filtered  HIMU 


Fig.  3.40.  Plot  of  F-test  significance  level  as  a  function  of  the  number  of 
singular  values  retained  for  the  filtered  HIMU  data  set.  Basically,  the  test 
determines  whether  additional  parameters  (singular  values]  make  a  significant 
contribution  to  the  model  fit  of  the  observed  data  values.  For  filtered  HIMU,  the 

optimal  p  value  [for  95%  significance]  is;  p  =  23.  Line  symbols:  o - o  =  filtered 

HIMU, - =  95%  significance  level. 
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Fig.  3.41.  Plot  of  F-test  significance  level  as  a  function  of  the  number  of 
singular  values  retained  for  the  DMM  data  set.  Basically,  the  test  determines 
whether  additional  parameters  [singular  values]  make  a  significant  contribution 
to  the  model  fit  of  the  observed  data  values.  For  DMM,  the  optimal  p  value  [for 

95%  significance]  is;  p  =  22.  Line  symbols:  o — o  =  DMM, - =  95% 

significance  level. 


Correlarion  Coefficient 
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Solutions  Minimizing  Coefficient  RMS  Error 


Degree 


Fig.  3.42.  Correlation  of  geophysics  coefficient  solutions,  that  minimize  the 

coefficient  rms  error,  with  the  actual  GEM-L2  coefficients.  Line  symbols: - 

=  geoid,  .  .  .  .  =  gravity,  -  .  -  .  =  gravity  gradient.  Confidence  levels  are 
determined  by  a  r-test  with  2/  degrees  of  freedom. 


Correlation  Coefficient 
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Solutions  Minimizing  Model  RMS  Error  -  Unfiltered 


Fig.  3.43.  Correlation  of  geophysics  coefficient  solutions,  that  minimize  the 
model  rms  error  for  the  unfiltered  data,  with  the  actual  GEM-L2  coefficients. 

Line  symbols; - =  geoid,  .  .  .  .  =  gravity,  -  .  -  .  =  gravity  gradient. 

Confidence  levels  are  determined  by  a  r-test  with  21  degrees  of  freedom. 
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Fig.  3.44.  Correlation  of  geophysics  coefficient  solutions,  that  minimize  the 
model  ims  error  for  the  filtered  data,  with  the  actual  GEM-L2  coefficients. 

Line  symbols; - =  geoid,  .  .  .  .  =  gravity,  -  .  -  .  =  gravity  gradient. 

Confidence  levels  are  determined  by  a  /-test  with  21  degrees  of  freedom. 
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CHAPTER  4 

SPHERICAL  HARMONIC  REPRESENTATION  OF  ISOTOPIC 
SIGNATURES:  THE  DELTA-FUNCTION  MODEL 

INTRODUCTION 

As  mentioned  in  Chapter  3,  the  delta-function  model  represents  the  OIB 
reservoir  as  a  series  of  point  sources,  each  feeding  a  separate  plume.  This  may 
seem  unphysical,  but  could  be  a  good  approximation  of  actual  conditions  if  the 
source  boundary  layer  is  not  continous,  but  patchy,  as  indicated  in  some  seismic 
studies  of  D"  (Lay  et  al.,  1990). 

Representing  the  geographic  features  as  delta-functions  [scaled  by  the 
corresponding  geoid  anomaly  or  mantle  component  percentage]  has  two 
advantages,  mathematically,  over  the  approximation  methods  used  in  Chapter  3. 
First,  the  spherical  harmonic  coefficients  can  be  found  easily  with  the 
simplification  from  integration  over  the  globe  to  summation  over  the  feature 
locations  allowed  by  the  delta-functions.  Second,  representing  the  OIB  reservoir 
as  a  known  function  removes  the  problem  of  aliasing;  the  values  of  the  spherical 
harmonic  coefficients  are  not  dependent  upon  the  truncation  point  of  the 
expansion  (they  are  dependent  upon  the  number  and  location  of  the  geographic 
features].  For  delta-functions,  which  have  energy  at  all  degrees,  the  expansions 
can  be  carried  out  to  infinity,  but  for  this  study,  will  only  be  carried  out  to 
degree  5,  for  comparison  witli  the  continuous  layer  model. 


THf:ORY 

As  before,  any  function /(G.tp)  can  be  expanded  in  spherical  harmonics: 


^G.(p) 


/,  / 

X  I 

/=<)  w  =  0 


/ 


A'  4  jr( /+//;)' 


I’/^VcosG)  1  A7'cos/?;(P  +  Bj'sin/utp] 


1  68 


Due  to  the  orthogonality  of  the  spherical  harmonics,  the  equations  for  the 


coefficients  are; 


aT^ 


I'A 


y(e,(p)  PTCcosB)  cos/>2(p  r^(cose) 


bT  = 


[  dA  y(e,(p)^/^^^SZ^rp^(cose)sinm(p£/(cose) 
j  J  V  4n(/+m)! 


For  the  delta-function  model,  the  function  being  expanded  is  a  series  of  delta- 


functions: 


/(e,(p)  =  ^,6(e-e„(j>-(p,) 


where  is  one  of  the  four  mantle  component  percentages  [or  the  value  of  the 

geoid  anomaly]  and  6(0-0i,(f)-(p,)  indicates  a  delta-function  at  the  particular 
location  (0/,(p/).  Mathematically,  the  delta-function  is  a  "spike"  of  infinite  height, 

infinitesimal  widlii  and  unit  area: 


j  d(pj  5(0-0(,(j)-(p,)  dB  =  1 

The  key  property  of  the  delta-function  is  that  the  integral  of  a  function  g(0,tp) 
times  a  delta-function  is  just  the  value  of  g  at  the  delta-function  location: 


d(p  g(0,(p)  5(0-0„(|)-(p,)  dS  =  g(0/,(p,) 


7'his  simplifies  the  coefficient  equations  from  integration  over  the  globe  to 
summation  over  the  geographic  feature  locations: 


a"'  =  X  p''’(cos0,)  coswtp, 

'  ^  V  47r(/+/n)! 
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N 

i=\ 


■  -rr  7,  P/CcosG/)  sinmcp, 
47c(/+m)! 


The  coefficient  equations  for  the  constructed  data  sets  of  geoid,  gravity 
and  gravity  gradient  anomalies  at  the  feature  locations  have  additional  factors. 
As  an  example,  for  gravity  the  equations  are: 


A 


m 

I 


GM(/+I)“  V  4)t(;+m)! 


P/”(cos0/)  cosmcp^ 


B 


m 

I 


y  /  (2/4-1  )(/-m)! 
GMil+\)^^  'y  4n(l+m)\ 


PT'CcosG/)  sinmcp; 


_R^ _ 

with  the  additional  factor  of  GM(\+\).  Geoid  and  gravity  gradient  additional 


factors  are  and  GM(l+lXl+2),  respectively. 


APPLICATION 

As  before,  the  constructed  geophysics  data  sets  are  used  as  a  control  to 
gauge  the  level  of  accuracy  expected  from  the  mantle  component  data  sets. 
Correlating  the  coefficients  from  these  data  sets  with  the  GEM-L2  coefficients 
(Fig.  4.1)  yields  good  agreement  for  all  three  at  degree  2.  Whereas  the 
continuous  layer  model  showed  a  fairly  consistent  pattern  of  decreasing 
correlation  from  the  geoid  coefficient  estimates  to  the  gravity  and  gravity 
gradient  estimates  (Fig.  3.43),  the  delta-function  model  shows  equal  correlation 
at  degree  2  and  a  switch  to  increasing  correlation  from  the  geoid  estimates  to  the 
gravity  and  gravity  gradient  estimates  at  degree  4.  Overall,  it  appears  that  the 
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delta-function  model  is  less  accurate  at  reproducing  the  coefficients  for  long 
wavelength  data  sets  Igeoid]  and  more  accurate  at  reproducing  the  coefficients 
for  the  short  wavelength  data  sets  [gravity  gradient]  than  the  continuous  layer 
model.  Both  models  are  consistent,  though,  in  showing  strong  correlation  for  all 
three  data  sets  at  degree  2,  implying  that  the  mantle  component  degree  2 
coefficients  are  also  viable.  In  addition,  the  mantle  component  data  sets  have 
even  more  high  degree  [short  wavelength]  energy  than  the  gravity  gradient  data 
set,  so  their  coefficients  are  probably  reasonably  accurate  out  to  degree  4. 

Since  each  of  the  different  geophysics  data  sets  approximate  the  GEM-L2 
coefficients  equally  well  at  degree  2,  it  appears  that  there  is  some  additional 
controlling  factor  affecting  the  estimates  of  the  degree  2  coefficients,  aside  from 
the  data  values  themselves.  The  location  of  the  features,  and  thus  the  delta- 
functions,  is  the  most  likely  candidate.  A  plot  of  the  constructed  degree  2 
"function"  for  the  delta-function  model  geoid  (Fig.  4.2)  shows  the  obvious 
relationship  between  the  two  main  clusters  of  oceanic  islands  and  the  two  highs  in 
the  geoid.  Since  the  continuous  layer  model  geophysics  coefficients  all  agreed 
well  w'ith  the  degree  2  GEM-L2  coefficients,  it  appears  that  the  location  effect 
merely  enhances  an  already  existing  correlation  and  is  not  solely  responsible  for 
the  correlation.  Presumably  the  same  is  true  of  any  degree  2  correlation  of 
delta-function  model  geochemistry  coefficients  with  the  GEM-L2  coefficients. 

Degree  2  "functions  '  for  the  mantle  component  percentages  are 
reconstructed,  as  before,  for  comparison  with  those  of  the  continuous  layer 
model  (Figs.  4. 3-4. 6).  The  contoured  values  of  the  delta-fuction  geoid  (Fig.  4.2) 
and  the  mantle  component  functions  are  large  enough  to  be  the  actual  geoid  and 
component  percentages,  instead  of  deviations  from  the  average  values,  as  for  the 
continuous  layer  model.  This  is  due  to  the  arbitrary  .scaling  that  comes  into  play 
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when  using  delta-functions.  A  delta-function  has  unit  area,  so  the  average  value 
of  a  delta-function  over  the  sphere  is: 

(6)  = - ^ =  ^ 

(Atp  sin0)  A0  471 

where  (A(psin0)  A0  is  a  sectional  area  on  the  sphere  (Fig.  4.7),  which  for  the 
whole  sphere  is  471.  If  there  is  only  one  delta-function  involved  in  the 
reconstruction,  the  contoured  values  will  be  off  by  a  factor  of  l/(47i).  Since 
there  are  36  features,  there  are  36  delta-functions  involved  in  the  reconstruction, 
so  the  contoured  values  are  off  by  a  factor  of  36/(47u)  =  2.86  or  ~3. 

Qualitatively,  the  four  reconstructed  mantle  component  degree  2  functions 
show  good  agreement  with  each  other.  All  four  have  two  highs:  one  over  central 
Africa  and  the  other  over  the  central  Pacific.  Slight  differences  include  the 
width  of  the  highs  [from  narrowest  to  widest  width:  HIMU,  EMI.  EMU  and 
DMM]  and  the  amount  of  displacement  [from  0°  to  15°]  of  the  highs  above  and 
below  the  equator  [from  least  to  most  displacement:  HIMU,  EMU,  EMI  and 
DMM].  With  respect  to  the  GEM-L2  degree  2  geoid  (Fig.  3.50),  all  of  the 
mantle  component  highs  are  shifted  longitudinally  to  the  east  by  varying  amounts 
[HIMU  -30°,  EMI  -30°,  EMU  -35°  and  DMM  -40°]. 

Degrees  2-3  functions  for  the  four  components  (Figs.  4.8-4. 1 1)  are 
constructed  for  comparison  with  the  geoid  (Fig.  5.32)  and  the  HIMU  continuous 
layer  model  reconstruction  (Fig.  5.31). 

SUMMARY 

Viewing  the  distribution  of  the  OIB  reservoir  as  a  series  of  point  sources 
that  can  be  represented  as  delta-functions  yields  the  following  results: 
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•  With  respect  to  the  behavior  of  geophysics  control  data  sets,  at  least  the 
degree  2  spherical  harmonic  coefficients  for  the  mantle  components 
can  be  estimated  with  confideiiee,  if  not  the  degrees  3  and  4  as  well. 

•  The  location  of  the  features,  and  thus  the  delta-functions,  biases  the 
calculated  degree  2  coefficients  due  to  the  correlation  between  the 
oceanic  island  locations  and  the  degree  2  geoid. 

•  Scaling  of  delta-function  models  reconstructed  over  the  globe  is 
dependent  upon  the  number  of  delta-functions  used  m  the 
approximation  [A^]  and  varies  as  A//(47t). 

•  Degree  2  HIMTJ,  EMI,  EMU  and  DMM  all  show  a  degree  2  geoid 
pattern  phase-shifted  30°-40°  to  the  east,  with  varying  widths  of  the 
highs  and  displacements  from  the  equator. 
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Delta-Function  Model  Geophysics  Correlation  with  GEM-L2 
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Fig.  4.1.  Correlation  of  the  delta-function  model  geophysics  coefficient  solutions 

with  the  actual  GEM-L2  coefficients.  Line  symbols: - =  geoid,  .  .  .  .  = 

gravity,  -  .  -  .  =  gravity  gradient.  Confidence  levels  are  determined  by  a  t-test 
with  21  degrees  of  freedom. 
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CHAPTER  5 

RESULTS  AND  DISCUSSION 

INTRODUCTION 

Geophysical  control  data  sets  are  used  to  judge  the  dependability  of 
spherical  harmonic  coefficient  solutions  for  the  mantle  end-member  components 
from  the  continuous  layer  and  the  delta-function  models.  A  careful  comparison 
of  the  two  models  can  further  enhance  or  reduce  the  significance  assigned  to  the 
various  solutions.  In  this  chapter,  the  two  models  are  compared  in  terms  of  their 
amplitude  spectra,  how  well  they  correlate  with  the  geoid,  how  they  are  affected 
by  nonuniform  feature  distribution  and  how  well  they  correlate  with  the 
Clayton-Comer  seismic  tomography  model.  The  significance  of  the  correlations 
with  the  geoid  and  the  seismic  tomography  model  is  discussed,  along  with 
suggestions  for  further  research. 

AMPLITUDE  SPECTRA 

Spectral  amplitude  plots  show  the  relative  power  at  each  degree  for  the 
different  mantle  component  expansions.  Following  Richards  and  Hager  (1988), 
the  root  mean  square  harmonic  coefficient  amplitude  at  each  degree  is  given  by: 


t  [(Ar!-+(Brp] 

71=0 _ 

{2/+I) 

where  is  the  variance  at  each  degree  for  a  given  set  of  harmonic  coefficients. 
Richards  and  Hager  (1988)  include  the  factor  of  1/(2/  +  1)  because  random  noise 
on  a  sphere  will  have  a  flat  spectrum  with  this  normalization.  On  plots  of 
versus  /,  low-degree  or  long-wavelength  effects  will  show  up  as  a  negative 


slope. 


1  96 


Anipliuidc  spectra  ot  the  calculated  geoid  coefficients  from  the  two 
models  agree  well  w'ith  the  negative  [long-w'avelength]  slope  of  the  actual  geoid 
coefficients  (I-igs.  5.1  and  5.2).  For  the  mantle  component  e.xpansions, 
amplitude  spectra  reveal  no  such  clear  cut  negative  slope  pattern  to  indicate 
dominant  long-wavelength  effects  (Figs.  5.3  and  5.4).  Instead,  the  spectra  appear 
"white  ",  with  energy  at  all  degrees,  and  no  decrease  in  the  energy  with 
increasing  degree.  In  addition,  HIMU  is  the  only  mantle  component  that  shows 
anv  consistency  in  behavior  betw'een  the  two  models.  Thus,  in  general,  the 
expansion  of  the  mantle  components  is  model  dependent. 

CORRFLATION  WITH  THE  GEOID 

Plotting  the  mantle  (  nniponent  percentages  point  by  point  against  the  full 
geoid  value  at  the  geographic  feature  locations  is  not  a  valid  way  to  compare  the 
mantle  component  signatures  with  the  geoid.  When  correlating  them  by  degree 
using  spherical  harmonic  coefficients,  it  is  apparent  that  the  mantle  components 
may  correlate  with  the  geoid  at  some  degrees  [wavelengths]  and  not  others.  In  a 
pointwise  comparison,  the  different  patterns  at  the  different  degrees  are 
obscured  as  they  are  added  together  to  produce  the  w'hole,  making  an  accurate 
comparison  impossible.  Pointwise  plots  done  with  the  current  data  show  no 
correlation  between  the  mantle  components  and  the  geoid  (Figs.  5. 5-5. 8). 

In  contrast,  correlating  the  geoid  coefficients  and  the  mantle  component 
coefficients  by  degree  reveals  a  good  corrrelation  (90%  significance  level  and 
higher]  at  degree  2  for  the  DUPAL  components  [EMI,  EMU  and  HIMU]  for  both 
models  (Figs.  5.9  and  5.10).  Note  that  positive  correlations  indicate  high 
concentrations  of  mantle  components  correlating  with  geoid  highs  and  vice  versa. 
HIMU  has  the  best  correlation  for  both  models,  showing  better  liian  95% 
significance  ;tt  degree  2  and  90%  significance  at  degree  3.  The  remaining  mantle 
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components  show  a  consistent  decreasing  correlation  from  EMU  to  EMI  to 
DMM  for  both  models. 

IMPLICATIONS  OF  NONUNIFORM  FEATURE  DISTRIBUTION 

Oceanic  island  distribution  is  not  uniform  about  the  globe.  As  indicated  in 
Chapter  4,  the  two  main  clusters  of  oceanic  islands  correspond  to  the  two  highs 
of  the  degree  2  geoid.  It  can  be  argued,  then,  that  any  correlation  between  the 
degree  2  mantle  component  expansions  and  the  degree  2  geoid  is  due  solely  to  the 
nonuniform  distribution  of  the  oceanic  islands  and  not  to  any  pattern  in  the 
geochemistry  values.  To  test  this,  the  percentages  of  the  HIMU  mantle 
component  at  the  36  geographic  features,  filtered  (continuous  layer  model  |  and 
unfiltered  (delta-function  model (,  are  randomly  assigned  to  different  feature 
locations  five  times.  HIMU  percentages  are  used  since  the  degree  2  HIMU,  for 
both  models,  correlates  best  with  the  degree  2  geoid.  The  five  randomly 
generated  data  sets  for  each  model  are  then  used  to  compute  new  coefficients  that 
can  be  compared  to  the  degree  2  geoid.  For  the  continuous  layer  model,  the 
number  of  singular  values  retained  for  the  new  data  sets  is  determined  by  the  F- 
test  at  the  95%  significance  level.  The  random  number  generator  used  for  this 
test  is  nonlinear,  but  repeatable,  since  it  starts  with  a  given  seed  that  is  updated 
for  successive  calls  in  a  predictable  manner.  This  means  that  for  a  given 
randomization,  the  filtered  and  unfiltered  HIMU  percentages  are  being 
randomized  in  the  same  way,  so  the  results  of  the  two  models  can  be  compared. 
Five  iterations  is  not  enough  to  quantify  the  effect  of  the  feature  distribution  on 
the  degree  2  correlation  for  the  two  models,  but  it  is  enough  to  indicate  if  it  has 
any  control  at  all. 

Concentrating  on  the  degree  2  coefficients,  three  of  the  randomizations 
that  result  in  strong  correlations  with  the  geoid  for  delta-function  model  (well 
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above  the  90%  conl'idence  level]  result  in  negligible  correlations  with  the  geoid 
for  the  continuous  layer  model  (Table  5.1).  Reconstructed  degree  2  functions  of 
the  randomized  data  sets  show  graphically  how  little  the  delta-function  model 
changes,  with  respect  to  the  continuous  layer  model,  when  the  geochemical 
signatures  of  the  features  are  mixed  up  (Figs.  5.1 1-5.20).  For  the  delta-function 
model,  this  indicates  that  the  values  of  coefficients  are  not  so  much  dependent 
upon  the  scaling  factors  multiplying  the  delta-functions  as  the  location  of  the 
delta-functions  themselves.  This  location  effect  makes  it  difficult  to  trust  strong 
correlations  of  the  delta-function  model  wdth  the  geoid  unless  there  is  additional 
confirmation  by  the  continuous  layer  model. 

CORRF.LATION  WITH  SEISMIC  TOMOGRAPHY 

Correlating  the  mantle  component  expansions  with  the  geoid  gives  an 
estimate  of  the  general  OIB  source  region  jie.  lower  mantle  versus  upper 
mantle],  but  is  incapable  of  resolving  a  more  precise  depth  range  for  the  source 
since  the  geoid  is  affected  by  mass  anomalies  at  all  depths  in  the  Earth.  A  way  to 
select  a  probable  depth  range  for  the  OIB  sourcejs]  is  to  compare  the  mantle 
component  expansions  to  seismic  tomography  models.  Seismic  tomography 
models  map  the  global  distribution  of  lateral  velocity  variations  in  the  mantle  at 
different  depths  based  upon  the  inversion  of  travel  time  anomaly  data  from 
seismic  waves  that  travel  through  the  Earth's  interior  (Hager  and  Clayton,  1989). 

In  this  study,  the  mantle  component  expansions  are  correlated  with  the 
Clayton-Comer  seismic  tomography  model,  discussed  in  Hager  and  Clayton 
(1989).  The  Clayton-Comer  model  inverts  for  slowness  [inverse  of  velocity] 
anomalies,  in  a  given  shell,  that  are  converted  to  velocity  anomalies  by 
multiplying  by  the  average  shell  velocity.  There  are  29  shells  in  the  model,  each 
100  km  thick,  spanning  the  entire  mantle  from  the  core-mantle  boundary  [CMB|, 
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at  a  depth  of'  2900  km,  to  the  surface.  Shells  23-29  [covering  the  uppper  mantle] 
are  not  used  in  this  analysis  since  coverage  in  the  top  700  km  of  the  mantle  is 
poor  because  of  the  near  vertical  seismic  ray  paths  in  this  region.  The  spherical 
harmonic  coefficients  of  the  remaining  22  shells  [covering  the  lower  mantle]  are 
averaged  together,  to  dampen  model  noise,  to  produce  5  layers:  2900-2500  km 
[layer  1  ],  2500-2100  km  flayer  2],  2100-1700  km  [layer  3],  1700-1200  km  [layer 
4]  and  1200-700  km  [layer  5]. 

3'he  gcoid  is  correlated  with  the  Clayton-Comcr  tomography  model  first 
(Fig.  5.211  to  serve  as  a  guide  for  interpreting  the  correlation  of  the  tomography 
model  with  the  mantle  component  expansions.  Note  that  a  negative  correlation 
indicates  geoid  highs  correlating  with  low  velocity  regions  [and  vice  versa]  and  a 
positive  correlation  indicates  geoid  highs  correlating  with  high  velocity  regions 
[and  vice  versa].  In  layers  1-3,  the  strong  negative  correlations  at  degrees  2  and 
3  confirm  that  long  wavelength  geoid  high.s  are  due  to  low  density  (warmer  and 
thus  slower  velocity]  mantle  upwellings.  This  long  wavelength  upwelling 
signature  is  also  present  in  the  upper  lower  mantle,  as  shown  by  the  strong 
negative  correlations  at  degrees  2,  3  and  4  for  layer  4  and  at  degree  2  for  layer 
5.  Of  interest  is  the  strong  positive  correlations  for  layers  4  and  5,  at  degree  5 
and  degrees  4  and  5,  respectively.  Bowin  (1991a)  indicates  the  correspondence 
of  the  degrees  4-10  geoid  highs  with  plate  convergence  zones.  He  believes  that 
the  mass  anomalies  responsible  for  the  highs  lie  in  the  lower  mantle,  beneath 
plate  convergence  zones,  below  the  teleseismically  downgoing  subducted  slabs. 
The  positive  correlations  in  layers  4  and  5  support  this  theory  and  imply  that 
subducted  slabs  extend  below  the  670  km  discontinuity. 

Correlation  of  the  mantle  component  expansions  with  the  Clayton-Comer 
tomography  layers  for  the  tw'O  models  yields  interesting  results  (Figs.  5.22- 
5.29).  Due  to  the  limitations  of  both  models  fie.  the  uncertainties  in  the 
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coeff  icient  estimates  for  the  continuous  layer  model  and  the  location  dependence 
in  the  delta-function  model],  it  is  more  likely  that  a  significant  correlation  is 
accurate  if  it  is  present  in  both  models.  With  this  in  mind,  the  interpretation  of 
the  correlation  results  will  be  based  upon  common  correlations  of  90% 
significance  [or  very  close  to  it]  or  higher  (Table  5.2). 

The  common  degree  2  correlations  with  layers  3-5  for  all  the  mantle 
components  are  indicative  of  large  scale  upwelling,  as  for  the  geoid.  Good 
degree  3  correlations  with  layer  1  points  to  a  deep  source  for  all  four 
components,  like  the  geoid  which  shows  a  much  stronger  correlation  at  degree  3 
with  layer  1  than  it  does  at  degree  2.  This  correlation  is  not  unexpected  for  the 
DIJPAL  components,  whose  correlation  with  the  degree  2  geoid  also  suggest  a 
deep  origin,  but  it  is  surprising  for  the  DMM  component.  There  are  two 
po'^sible  solutions  for  the  dilemma  po.sed  by  the  supposedly  upper  mantle  DMM 
component  correlating  with  deep  mantle  tomography.  First,  it  is  possible  that 
the  DMM  component  expansion  does  correlate  better  with  upper  mantle 
tomography,  which  is,  unfortunately,  not  available  for  the  Clayton-Comer 
model.  Second,  it  is  possible  that  the  DMM  component  is  representative  of  both 
the  upper  and  lower  mantle  composition.  Han  (1991)  shows  that  all  the  hotspots 
that  have  elongated  isotopic  arrays  indicate  mixing  between  one  of  the  DUPAL 
components  and  something  that  is  not  a  MORB  composition.  Since  3/4  of  a 
plume's  ascent  is  spent  in  the  lower  mantle,  the  composition  of  the  DMM 
component  may  be  largely  controlled  by  lower  mantle  entrainment  (Hart,  1991). 

Another  interesting  correlation  common  to  both  models  is  the  positive 
correlation  at  degree  5  for  EMII  in  layer  5.  With  respect  to  Bowen's  model 
(1991a)  this  indicates  a  correlation  between  the  EMU  component  and  subducted 
slabs.  This  finding  agrees  with  the  geochemical  evidence  suggesting  the  EMU 
component  is  derived  from  recycling  of  subducted  sediments  (Hart,  1988). 
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DISCUSSION 

As  indicated  in  Chapter  3,  the  average  value  of  the  geoid  anomaly  at  the 
36  feature  locations  is  13.7  m,  not  zero  as  it  should  be  if  the  features  were 
located  randomly  with  respect  to  the  geoid.  This  is  a  simple  indication  that  the 
feature  locations  (hotspots]  correlate  with  geoid  highs.  Naturally,  then,  the  bulk 
chemical  signatures  unique  to  oceanic  island  basalts  should  also  correlate  with 
geoid  highs.  What  is  significant  is  that  the  expansions  of  all  three  DUPAL 
mantle  end-member  components  (EMI,  EMU  and  HIMU],  that  comprise  3/4  of 
the  bulk  chemical  signature,  individually  correlate  with  geoid  highs.  More 
importantly,  the  DUPAL  components  correlate  with  the  degree  2  geoid  highs, 
indicating  a  deep  origin  for  the  components  since  the  degrees  2-3  geoid  field  is 
inferred  to  result  from  topography  at  the  core-mantle  boundary  (Bowen,  1991a). 

It  can  be  argued  that  the  correlation  of  the  DUPAL  components  with  the 
degree  2  geoid  is  not  an  indication  of  geochemical  patterns  within  the  earth,  but  a 
direct  result  of  the  nonuniform  distribution  of  the  oceanic  islands,  whose  two 
largest  population  densities  correspond  to  the  degree  2  geoid  highs. 
Randomization  tests  indicate,  however,  that  while  this  nonuniform  distribution 
does  play  a  role  in  solutions  for  the  delta-function  model,  it  is  not  the  controlling 
factor  for  continuous  layer  model  solutions.  Though  the  continuous  layer  model 
solutions  are  hindered  by  the  limited  number  and  coverage  of  the  oceanic  islands 
and  the  delta-function  model  solutions  are  biased  by  the  oceanic  island  locations, 
continual  comparisons  of  the  two  models  can  be  used  to  judge  the  accuracy  of  the 
solutions  (in  addition  to  judging  accuracy  using  geophysical  control  sets]. 
Essentially,  where  both  models  agree,  the  solutions  are  more  likely  to  be 


accurate. 
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riie  total  gcoid  field  is  due  to  the  contribution  of  different  mass  anomalies 
at  difierent  depths  throughout  the  Earth,  so  it  can  be  difficult  to  directly 
ascertain  a  source  depth  by  comparing  geochemical  quantities  with  the  geoid. 
Seismic  tomography  models  allow  the  correlation  of  geochemical  quantities  with 
seismic  velocity  anomalies  at  different  depths  and  serve  as  an  independent  check 
on  the  general  source  locations  indicated  by  correlation  with  geoid  anomalies. 
Correlating  the  mantle  end-member  components  from  both  models  with  the 
Clavton-Comer  seismic  tomography  model  suggests  a  source  depth  range  of 
25()()-29()()  km  |just  above  the  core-mantle  boundary]  for  the  DIJPAL 
components,  due  to  the  strong  negative  degree  3  correlations  at  this  depth.  In 
addition,  a  strong  positive  degree  5  correlation  in  the  depth  range  of  700- 1200 
km  is  an  indication  that  the  EMll  component  is  related  to  subduction,  as 
previously  suggested  using  geochemical  evidence  (Hart,  1988).  Similarly,  the 
geoid  shows  a  strong  positive  correlation  with  the  Clayton-Comer  model  at 
degrees  4  and  5  in  the  depth  range  700-1200  km  and  at  degree  5  in  the  depth 
ranges  of  12(X)-1700  km.  These  subduction  related  patterns  in  the  upper  lower 
mantle  indicate  that  subducted  slabs  extend  beyond  the  670  km  seismic 
discotitinuity  and  thus  are  supporting  evidence  for  whole  mantle  convection 
b'urther  comparisons  need  to  be  made  between  the  mantle  component 
expansions  and  other  seismic  tomography  models.  It  is  especially  important  to 
ccM'ipare  the  mantle  components  to  a  high  resolution  upper  mantle  tomography 
model,  since  the  amplitude  spectra  for  the  components  indicate  power  at  high 
decrees  which  will  become  dominant  at  shallow  depths  in  the  mantle.  Such  a 
comparision  could  clarify  the  nature  of  the  DMM  component,  which  correlates 
well  with  the  degree  3  deep  mantle  layer  of  the  Clayton-Comer  model,  and  could 
further  explore  the  relationship  betw-een  the  EMU  component  and  subduction. 
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SUMMARY 

A  comparison  of  the  two  models  used  to  expand  the  mantle  components  in 
spherical  harmonics  yields  the  following  results: 

•  Mantle  end-member  component  amplitude  spectra,  for  the  continuous 
layer  model  and  the  delta-function  model,  show  power  at  all  degrees, 
with  no  one  degree  dominating. 

•  The  DUP.AL  components  (HMI,  HMII  and  HIMLJJ  for  both  models 
correlate  well  with  the  geoid  at  degree  2,  indicating  a  deep  origin. 

•  Dcltii-function  model  solutions  arc,  to  some  extent,  controlled  by  the 
nonuniform  feature  distribution,  while  the  continuous  layer  model 
solutions  arc  not. 

•  The  DUPAL  and  DMM  components,  for  both  models,  correlate  well 
[negatively!  at  degree  3  with  the  velocity  anomalies  of  the  Claytott- 
Comer  seismic  tomography  model  in  the  2500-2900  km  depth  range 
(immediately  above  the  core-mantle  boundary]. 

•  The  13M11  component,  for  both  models,  correlates  well  [positively]  at 
degree  5  with  the  velocity  anomalies  of  the  Clayton-Comer  seismic 
tomography  model  in  the  700-1200  km  depth  range,  indicating  a 
subduction  related  origin. 

•  .Subduction  related  positive  correlations  for  the  geoid  and  the  EMU 
component  with  the  Clayton-Comer  model  in  the  upper  lower  mantle 
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1  /(X)  17(K)  kin  I  indicaic  that  sulxluctcd  slabs  exiciul  bolnvv  the  670  kin 
seismic  discontinuity,  supporting  a  wliolc-inantle  convection  motlel. 
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Table  5.1.  Siiiniiiar\'  of  correlalion  coefficients  between  the  GliM-L2 
coefficicnis  and  coefficients  calcul:  ted  from  five  randomly  generated  data  sets 
for  the  continuous  layer  model  |filt(‘red  HIMUJ  and  the  delta-function  model 
llll.Ml.'l,  along  with  the  actual  corrcladons  of  the  filtered  HIMU  and  HlMl  ’  data 

sets. 


Data  Set 

Degree  2 

Degree  3 

Degree  4 

Degree  5 

(ontinuous  I. aver 

Model 

filtered  IIIMII 

0.752 

0.502 

-0.1 12 

-0.358 

raiuh'in  1 

0.753 

0.446 

-0.639 

-0.157 

raiKhun  ? 

0.560 

0.196 

0.467 

-0.210 

random 

-0. 1 29 

0.432 

-0.303 

-0.06b 

random  4 

0.225 

0.448 

0.386 

-0.130 

randi'm  5 

-0.166 

0.718 

-0.285 

-0.230 

Delta-1- Unction  .Model 

UlMi: 

0.850 

0.491 

0.063 

-0.505 

random  1 

0.726 

0.332 

0.036 

-0.416 

random  2 

0.622 

0.404 

0.107 

-0.361 

random  3 

0.873 

0.477 

0.029 

-0.320 

random  -i 

0.893 

0.407 

0.415 

-0.385 

random  5 

0.761 

0.383 

0.107 

-0.286 
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rable  5.2.  Summary  of  correlations  of  90%  significance  [or  very  close  to  it]  or 
higher  for  the  continuous  layer  model  and  the  delta-function  model  when 
correlated  with  five  averaged  layers  in  the  Clayton-Comer  tomography  model. 
[A  or  next  to  the  component  name  indicates  a  positive  or  negative 
correlatioti,  respectively. j 


Degree  2 

Degree  3  Degree  4  Degree  5 

l.aycr  5 

-EMI 

-EMII 

-fllMU' 

+EMII 

I. aver  4 

-EMI 

-EMll' 

-dmm' 

Layer  3 

-EMI 

-EMIl' 

Layer  2 

-dmm2  -emu' 

I. aver  1 

-EMi' 

-EMI 

-EMll'  -EMIl' 

-HIMU 

-DMM  -DMM 

'The  continuous  layer  model  correlation 
^The  delta-function  model  correlation  is 

is  slightly  less  than  90%  significant, 
slightly  less  than  90%  significant. 

Degree 


Fig.  5.1.  Amplitude  spectra  for  the  continuous  layer  model  coefficients  of  the 

constructed  geoid  data  set,  as  compared  to  the  actual  geoid.  Line  symbols: - 

-  constructed  geoid, - =  GEM-L2  geoid. 


RMS  Amplitude 


xiu  ’ 

2.5, - ---T- 


Delta-Function  Model  Geoid 


Degree 


Fig.  5.2.  Amplitude  spectra  for  the  delta-function  model  coefficients  of  the 
constructed  geoid  data  set,  as  compared  to  the  actual  geoid.  Line  symbols: 
=  constructed  tteoid, - =  GEM-L2  geoid. 


RMS  Ampliuide 
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Fig.  5.3.  Amplitude  spectra  for  the  continuous  layer  model  coefficients  of  the 

mantle  component  data  sets.  Line  symbols: - =  filtered  EMI, - = 

LMll,  •  =  filtered  HH'/TU,  -  •  •  =  DMM. 
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Delta-Function  Model  Mantle  Components 


Fig.  5.4.  Amplitude  spectra  for  the  delta-function  model  coefficients  of  the 

mantle  component  data  sets.  Line  symbols: - =  EMI, - =  EMU,  •  •  •  •  = 

HIMU,  -  -  •  =  DMM. 
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%EMI 


Fig.  5.5.  Pointwise  comparison,  at  each  geographic  feature,  of  the  full  geoid 
anomaly  fin  meters]  with  the  EMI  component  percentage.  This  plot  gives  the 
impression  that  there  is  no  correlation. 


212 


%EM|| 


Fig.  5.6.  Poiniwise  comparison,  at  each  geographic  feature,  of  the  full  geoid 
anomaly  [in  meters]  with  the  EMU  component  percentage.  This  plot  gives  the 
impression  that  there  is  no  correlation. 
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%HIMU 


Fig.  5.7.  Pointwise  comparison,  at  each  geographic  feature,  of  the  full  geoid 
anomaly  [in  meters]  with  the  HIMU  component  percentage.  This  plot  gives  the 
impression  that  there  is  no  correlation. 


Continuous  Layer  Model  Mantle  Components  Correlation  with  GEM-L2 


Degree 


Fig.  5.9.  Correlation  of  the  continuous  layer  model  mantle  component 

coefficient  solutions  with  the  GEM-L2  geoid  coefficients.  Line  symbols: - : 

filtered  EMI, - =  EMH,  ■  •  •  •  =  filtered  HIML,  -  •  -  •  =  DMM.  Confidence 

levels  are  determined  by  a  r-test  with  21  degrees  of  freedom. 


Correlation  Coefficient 


Fig.  5.10.  Correlation  of  the  delta-function  model  mantle  component  coefficient 

solutions  with  the  GEM-L2  geoid  coefficients.  Line  symbols: - =  EMI,  — 

-  -  EMU,  ■  •  •  ■  =  HIMU,  -  ■  -  •  =  DMM.  Confidence  levels  are  determined  by  a 
/-test  with  21  decrees  of  freedom. 
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Coixelation  of  Geoid  with  Clayton-Comer  Model 


Degree 


Fig.  5.21.  Correlation  of  the  GEM-L2  geoid  coefficients  with  the  five  layers  of 

the  Clayton-Comer  seismic  tomography  model.  Line  symbols: - =  layer  1 

[2500-2900  ;:m], - =  layer  2  (2100-2500  km],  •  ■  ■  •  =  layer  3  (1700-2100 

km],  -  ■  -  •  =  layer  4  (1200-1700  km],  o - o  =  layer  5  [700-1200  km]. 

Confidence  levels  are  determined  by  a  f-test  with  21  degrees  of  freedom. 
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Continuous  Layer  Model  EMI  Correlated  with  Clayton-Comer  Model 


Fig.  5.22.  Correlation  of  the  continuous  layer  model  filtered  EMI  coefficients 
with  the  five  layers  of  the  Clayton-Comer  seismic  tomography  model.  Line 

symbols: - =  layer  1  [2500-2900  km), - =  layer  2  [2100-2500  km],  ■  ■  • 

•  =  layer  3  [1700-2100  km],  -  •  -  •  =  layer  4  [1200-1700  km],  o - o  =  layer  5 

[700-1200  km].  Confidence  levels  are  determined  by  a  Mest  with  21  degrees  of 
freedom. 


Continuous  Layer  Model  EMII  Correlated  with  Clayton-Comer  Model 


Fig.  5.23.  Correlation  of  the  continuous  layer  model  EMII  coefficients  with  the 
five  layers  of  the  Clayton-Comer  seismic  tomography  model.  Line  symbols: 

- =  layer  1  12500-2900  km], - =  layer  2  [2100-2500  km], - =  layer 

3  [1700-2100  km],  -  •  -  •  =  layer  4  [1200-1700  km],  o — o  =  layer  5  [700-1200 
km].  Confidence  levels  are  determined  by  a  /-test  with  2/  degrees  of  freedom. 
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Fig.  5.24.  Correlation  of  the  continuous  layer  model  filtered  HIMU  coefficients 
with  the  five  layers  of  the  Clayton-Comer  seismic  tomography  model.  Line 

symbols: - =  layer  1  [2500-2900  km], - =  layer  2  [2100-2500  km],  ■  ■  • 

•  =  layer  3  [1700-2100  km],  -  ■  -  •  =  layer  4  [1200-1700  km],  o - o  =  layer  5 

[700-1200  km].  Confidence  levels  are  determined  by  a  t-test  with  21  degrees  of 
freedom. 
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Delta-Function  Model  EMI  Correlated  with  Clayton-Comer  Model 
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Fig.  5.26.  Correlation  of  the  delta-function  model  EMI  coefficients  with  the  five 

layers  of  the  Clayton-Comer  seismic  tomography  model.  Line  symbols: - = 

layer  I  [2500-2900  km], - =  layer  2  [2100-2500  km],  •  •  •  •  =  layer  3  [1700- 

2100  km],  -  •  -  ■  =  layer  4  [1200-1700  km],  o - o  =  layer  5  [700-1200  km]. 

Confidence  levels  are  determined  by  a  Mest  with  21  degrees  of  freedom. 
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Delta-Function  Model  EMU  Correlated  with  Clayton-Comer  Model 
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Fig.  5.27.  Correlation  of  the  delta-function  model  EMII  :oefficients  with  the 
five  layers  of  the  Clayton-Comer  seismic  tomography  model.  Line  symbols; 

- layer  1  12500-2900  km], - =  layer  2  [2100-2500  km],  ■  •  •  •  =  layer 

3  [1700-2100  km],  -  •  -  •  =  layer  4  [1200-1700  km],  o — o  =  layer  5  [700-1200 
km].  Confidence  levels  are  determined  by  a  t-test  with  21  degrees  of  freedom. 
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Delta-Function  Model  HIMU  Correlated  with  Clayton-Comer  Model 


Fig.  5.28.  Correlation  of  the  delta-function  model  HIMU  coefficients  with  the 
five  layers  of  the  Clayton-Comer  seismic  tomography  model.  Line  symbols; 

- layer  1  12500-2900  km], - =  layer  2  [2100-2500  km],  ■  ■  •  •  =  layei 

3  [1700-2100  km],  -  •  -  ■  =  layer  4  [1200-1700  km],  o — o  =  layer  5  [700-1200 
km].  Confidence  levels  are  determined  by  a  Mest  with  21  degrees  of  freedom. 
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Fig.  5.29.  Correlation  of  the  delta-function  model  DMM  coefficients  with  the 
five  layers  of  the  Clayton-Comer  seismic  tomography  model.  Line  symbols: 

- =  laver  1  [2500-2900  km], - =  layer  2  [2100-2500  km],  •  •  •  ■  =  layc 

3  [1700-2100  km],  -  ■  -  =  layer  4  [1200-1700  km],  o — o  =  layer  5  [700-1200 
km].  Confidence  levels  are  determined  by  a  /-test  with  21  degrees  of  freedom. 
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